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ABSTRACT 


The sliding mode control is known as a robust control that is able to work under 
the uncertainties of modeling error and the environmental disturbances. The objective of 
this research is to design the simple control (sliding control) algorithms for a single-link 
flexible arm and to study the robustness due to varying payload. A general form of 
physical plant in state space is formulated. To achieve a continuous control, a time- 
varying boundary layer was introduced into the control system neighboring the sliding 
surface. The computer simulation program was coded in MATLAB. A low-cost IBM-AT 
micro-computer was utilized to implement the sliding control on the flexible arm system. 


A . i m» i i S fa ; 

■ ’>:<»*l 

r** T*t' 

.v (jut »a 

! Jtiruist 1 <jn 


i 


»/ - 
■■ r I Do. 

* * » ft c * . it 




TABLE OF CONTENTS 


l. INTRODUCTION. 1 

A. MOTIVATION . 1 

B. LITERATURE REVIEW . 2 

C. OBJECTIVE. 3 

II. PLANT AND ITS MATHEMATICAL MODELS. 4 

A. PHYSICAL PLANT . 4 

B. MATHEMATICAL MODEL OF THE PLANT. 5 

1. AN ERLS DYNAMICS MODEL OF FLEXIBLE ARM . . 5 

2. EQUATIONS OF MOTION OF FLEXIBLE ARM. 7 

3. ELECTROHYDRAULIC ACTUATION. 11 

C. STATE SPACE REPRESENTATION. 12 

m. CONTROLLER DESIGN (SLIDING MODE CONTROL) . 15 

A. INTRODUCTION. 15 

B. SYSTEM EQUATIONS AND CONFINED UNCERTAINTIES . . 17 

C. SLIDING SURFACES. 20 

D. SLIDING CONTROL WITH A FIRST-ORDER SLIDING 

CONDITION. 21 

E. STRAIGHT SLIDING CONTROL WITH FIRST-ORDER PLUS 

INTEGRAL SLIDING CONDITION. 30 

F. VERSATILE SLIDING CONTROL WITH SECOND-ORDER 


IV 





















SLIDING CONDITION 


36 


IV. RESULTS.44 

A. RESULTS ON SLIDING CONTROL WITH FIRST-ORDER 

SLIDING CONDITION PLUS INTEGRAL ERROR . 45 

B. RESULTS ON STRAIGHT SLIDING CONTROL. 47 

C. RESULTS ON VERSATILE SLIDING CONTROL. 50 

V. CONCLUSIONS AND RECOMMENDATIONS . 84 

A. CONCLUSIONS. 84 

B. RECOMMENDATIONS. 86 

APPENDIX A : SIMULATION PROGRAM . 87 

APPENDIX B : IMPLEMENTATION PROGRAM. 106 

LIST OF REFERENCES . 116 

INITIAL DISTRIBUTION LIST. 117 


v 
















LIST OF FIGURES 


PAGE 


1 A Synthesis View of Sliding Mode Control. 3 

2 Single-Link Flexible Arm. 5 

3 Generalized Coordinate. 6 

4 Block Diagram of Sliding Mode Control. 17 

5 A Detailed Synthesis View of Sliding Mode Control with First- 

Order Sliding Condition. 29 

6 A Detailed Block Diagram of Sliding Mode Control with 

First-Order Sliding Condition . 29 


7 A Detailed Synthesis View of the Straight Sliding Control . 37 

8 A Detailed Block Diagram of the Straight Sliding Control . 37 

9 A Detailed Synthesis View of the Versatile Sliding Control 43 

10 A Detailed Block Diagrm of the Versatile Sliding Control . 43 


11 (a)-11(d) The Time History of System Parameters 0.00 kg. 52 

11(e)-11(h) The Time History of System Parameters 0.85 kg. 53 


12(a)-12(d) The Simulation Performance of SMC with First-Order Sliding 

Condition (1^ = 10 (rad/sec), a t = 0.001 (sec), Payload = 

0.00 kg) . 54 

13(a)-13(d) The Simulation Performance of SMC with First-Order Sliding 

Condition (X^ = 500 rad/sec, a t = 0.001 sec, Payload = 0.00 


vi 














14(a)-14(d) 


kg). 55 

The Simulation Performance of SMC with First-Order Sliding 
Condition (1^ = 500 (rad/sec), At = 0.0005 (sec), Payload = 


0.00 kg) . 56 

14(e)-14(h) The Experiment Performance of SMC with First-Order Sliding 

Condition (X^ = 500 (rad/sec), Payload = 0.00 kg) ... 57 
15(a)-15(d) The Simulation Performance of SMC with First-Order Sliding 

Condition (X^ = 500 (rad/sec), At = 0.0005 (sec), Payload = 

0.85 kg) . 58 

15(e)-15(h) The Experiment Performance of SMC with First-Order Sliding 

Condition (X^ = 500 (rad/sec), Payload = 0.85 kg) ... 59 

16(a)-16(d) The Simulation of Tracking Performance of SMC with First-Order 

Sliding Condition (X,^ = 1700 (rad/sec), At = 0.0005 (sec). 

Payload = 0.00 kg. 60 

16(e)-16(h) The Experiment of Tracking Performance of SMC with First- 

Order Sliding Condition (X^ = 1700 (rad/sec), Payload = 0.00 

kg). 61 

17(a)-17(d) The Simulation of Tracking Performance of SMC with First-Order 

Sliding Condition (X^ = 1700 (rad/sec), At = 0.0005 (sec) 

Payload = 0.85 kg) . 62 

The Experiment of Tracking Performance of SMC with First- 

vii 


17(e)-17(f) 








Order Sliding Condition (X,^ = 1700 (rad/sec) Payload = 0.85 


kg). 63 

18(a)-18(d) The Simulation Performance of Straight Sliding Control 

(«„ = 10 (rad/sec), C t = 1, At = 0.0002 (sec), Payload = 0,00 

kg). 64 

19(a)-19(d) The Simulation Performance of Straight Sliding Control 

(«„ = 10 (rad/sec), C L = 10, Af = 0.0002 (sec), Payload = 0.00 

kg). 65 

19(e)- 19(h) The Experiment Performance of Straight Sliding Control 

(«. = 10 (rad/sec), ( L = 10, Payload = 0.00 kg). 66 

20(a)-20(d) The Simulation Performance of Straight Sliding Control 

(«„ = 1 (rad/sec), C t = 10, a t = 0.0002 (sec), Payload = 0.00 

kg). 67 

21(a)-21(d) The Simulation Performance of Straight Sliding Control 

(w, = 1 (rad/sec), = 100, At = 0.0002 (sec), Payload = 0.00 

kg). 68 

21(e)-21(h) The Experiment Performance of Straight Sliding Control 

(u„ = 1 (rad/sec), = 100, Payload = 0.00 kg). 69 

22(a)-22(d) The Simulation Performance of Straight Sliding Control 

(«, = 1 (rad/sec), C L = 100, At = 0.0002 (sec), Payload = 0.85 
kg). 70 

viii 











22(e)-22(f) 

23(a)-23(d) 

23(e)-23(f) 

24(a)-24(d) 

24(e)-24(f) 

25(a)-25(d) 

25(e)-25(h) 

26(a)-26(d) 

26(e)-26(h) 


The Experiment Performance of Straight Sliding Control 

(«„ = 1 (rad/sec), C t = 100, Payload = 0.85 kg). 71 

The Tracking Performance in Simulation of Straight Sliding 
Control («, = 1 (rad/sec), C L = 100, a t — 0.0002 (sec), 

Payload = 0.00 kg) . 72 

The Tracking Performance in Experiment of Straight Sliding 
Control = 1 (rad/sec), C L - 100, Payload = 0.00 kg) 73 
The Tracking Performance in Simulation of Straight Sliding 
Control (u, = 1 (rad/sec), C £ = 100, Af = 0.0002 (sec), 

Payload = 0.85 kg) . 74 

The Tracking Performance in experiment of Straight Sliding 
Control (u, = 1 (rad/sec), ( £ = 100, Payload = 0.85 kg) 75 
The Simulation Performance of Versatile Sliding Control 
(w„ = 1 (rad/sec), z - 1 C L = 10, At - 0.0002 (sec), Payload 

= 0.00 kg). 76 

The Experiment Performance of Versatile Sliding Control 
(<*>„ = 1 (rad/sec), * = 1 C £ = 10, Payload = 0.00 kg) . 77 

The Simulation Performance of Versatile Sliding Control 
(w,, = 1 (rad/sec), z - 10, ( L = 10, a t = 0.0002 (sec), Payload 

= 0.00 kg). 78 

The Experiment Performance of Versatile Sliding Control 


ix 








(«, = 1 (rad/sec), z = 10, C L = 10, Ar = 0.0002 (sec), Payload 


27(a)-27(d) 


27(e)-27(h) 


28(a)-28(d) 


28(e)-28(h) 


= 0.00 kg). 79 

The Simulation Performance of Versatile Sliding Control 
(u„ = 1 (rad/sec), z = 10, C L = 10, a* = 0.0002 (sec), Payload 
= 0.85 kg). 80 

The Experiment Performance of Versatile Sliding Control 
(«„ = 1 (rad/sec), z = 10, C L — 10, Payload = 0.85 kg) 81 
The Tracking Performance in Simulation of Versatile Sliding 
Control (u, = 1 (rad/sec), z = 10, C L = 10, a t = 0.0002 (sec), 

Payload = 0.00 kg) . 82 

The Tracking Performance in Experiment of Versatile Sliding 

Control (w, = 1 (rad/sec), z = 10, C x = 10, Payload = 0.00 
kg). 83 


x 








ACKNOWLEDGMENT 


First, I would like to give my thanks to Jesus Christ, my Lord, with his 
strengthening, I can have my thesis completed, and then, I wish to express my deepest 
gratitude to my advisor, Professor Liang-Wey Chang, for his unfailing guidance and full 
support for this research. Last but not least, I wish to thank my wife, Roselinda, and 
my son, Fourier, for their patience and help during my studies. 



I. INTRODUCTION 


A. MOTIVATION 

The control of flexible arms has been an active research area and a challenge to 
researchers. The motion control of flexible arms concerns the tip position and the control 
algorithm should be able to deal with the large motion and the small motion due to vibrations. 
It is difficult to have a precise arm model and an error exists between the model and the plant, 
that is modeling error. The modeling error for flexible arm may contain changing payload, high- 
frequency unmodeled dynamics, deviation on load position in the end-effector, and 
environmental disturbances. 

As the operation of control system is concerned, the computation of the control law plays 
an important role. While the time-delaying may worsen the control system, the control structure 
has to be as simple as possible to reduce the computation time. Since the models of flexible arms 
are complicated, simplified models will be utilized and controllers will be designed accordingly. 

Therefore, the challenge that we encounter in this research is to design simple control 
algorithms such that the modeling errors due to changing payload are compensated and the on¬ 
line operation will also be achieved on the single-link flexible arm. An IBM-AT is chosen for 
the low-cost implementation in this research. 

B. LITERATURE REVIEW 

The single-link flexible arm model of using Equivalent Rigid Link System (ERLS) was 
first derived by Chang [Ref. 1], The ERLS described the motion in large motion and small 
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motion. A Sequential Integration Method [Ref. 2] was also derived to facilitate an effective 
integration routine. Petroka [Ref. 3] experimentally validated the ERLS dynamic model of 
single-link flexible arm which was built and driven by an electrohydraulic actuator. Gannon 
[Ref. 4] upgraded the model by using the natural-mode shape function. Park [Ref. 3] designed 
and simulated a closed-loop non-robust controller for the arm. Kirkland [Ref. 5] redefined and 
implemented the controller on an IBM-AT computer. A strain gage and a potentiometer were 
used to determine the tip position of the arm. 

As the robust control algorithm is concerned, the Sliding Mode Control (SMC) [Ref. 6-9] 
has been known as a robust control that is insensitive to modeling error and disturbances. Figure 
1 illustrates a synthesis view of the sliding mode control. The external disturbances and 
parameters variation will first be filtered by the S dynamics to produce S and the generated S 
will be fed into error dynamics to perform further filtering. Thus, the S dynamics play very 
important roles to provide robustness to the control system. Once the S dynamics reject all the 
unwanted signals, the error dynamics will present the system behavior. Because of the 
undesirable high speed switched control on the sliding mode, a boundary layer thickness was 
introduced to the control system and a smoothed control was designed. Fan [Ref. 10] simulated 
and implemented the sliding mode control on the single link flexible arm. The robustness of the 
control was proved despite of a simplified model. Straight sliding control [Ref. 11] and versatile 
sliding control [Ref. 12] utilizing the idea of filtering provide more tuning capability on the S 
dynamics such that the unwanted signals will be filtered. In this research, these sliding control 
algorithms will be applied to the single link flexible arm to study the control system dynamics. 
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C. OBJECTIVE 


The objective of this research is to utilize simple control structures (sliding mode 
controls) to perform the motion control of a single-link flexible arm. The robustness due to 
varying payload will be studied. The control design will be performed in a state-space form 
using matrix-norm techniques. An IBM-AT will be chosen for the implementation. 
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H. PLANT AND ITS MATHEMATICAL MODELS 


A. PHYSICAL PLANT 

An experimental arm is shown in Fig. 2, which is driven by an electrohydraulic actuator. 
The motion of the arm is limited to a vertical plane. The flexible arm can bend freely in the 
vertical plane, but it is stiff in torsion and horizontal bending. The arm consists of two parallel 
steel flat bars welded at the base and directly clamped to the hydraulic actuator. Torsional 
stiffness is obtained by connecting the two steel bars, with thin steel strips, to seven transverse 
steel bridges. Table 1 shows geometric and mass properties of the flexible arm. 

Table 1 


Arm Length 

0.9985 m 

Arm Mass 

4.8565 kg 

Transverse Rigidity 

81.3 N m 2 

Arm Cross-Sectional Area 

6.178x10"* m 2 

Density 

7861.05 kg/m 3 


A potentiometer is used to measure the actuator (large motion) signal. A two-arm bridge 
strain gage attached to the center of the arm is used to calculate the tip deflection (small motion) 
of the arm. Data acquisition was performed using a high speed Data Translation interface board 














DT 2821-F-8DI, which was installed in a micro-computer (standard IBM-PC AT). The support 
software (AT-LAB) allowed direct manipulation of the data acquisition board through the use 
of provided subroutines which are compatible with FORTRAN. 



Figure 2 A Single-Link Flexible Manipulator System 


B. MATHEMATICAL MODEL OF THE PLANT 

In this section, a mathematical model of the plant is given, which includes an ERLS 
(Equivalent Rigid Link System) dynamic model of the arm and actuator dynamics. 

I. AN ERLS DYNAMIC MODEL OF FLEXIBLE ARM 

The ERLS is defined as a hypothetical system which produces the large motion 
and whose kinematics are equivalent to a rigid-link system. An ERLS of a planar manipulator 
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with single link is shown schematically with the dash line in Figure 3, where geometric center 
line of link was drawn. The solid lines stand for the deformed state of the arm. The ERLS of 
a flexible arm describes the large motion of the arm, and then the small motion arising from the 
structure flexibility can be superimposed on the ERLS [Ref. 11 

The ERLS dynamic model of the arm was developed by means of Lagrange’s 
formulation, the Finite Element Method, and the ERLS kinematics. To apply the Lagrangian 
dynamics to the flexible arm, the generalized coordinate (Figure 3) is chosen to describe the 
large motion by joint variable 0. The Finite Element Method is utilized to discretize the 
displacement such that the small motion is represented in terms of nodal displacement u , where u 
is the tip deflection. In this study, the natural mode shape functions of a beam are used to 
represent the flexural motion of the flexible arm, and the arm is modeled as a continuous Euler- 
Bemoulli cantilever beam, neglecting shear deformation and rotary effects. [Ref. 4] 
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2. 


EQUATIONS OF MOTION OF FLEXIBLE ARM 


Lagrange’s equations for the flexible manipulator are 


d.dKE s dKE dPE 
—(-)-+-=GF. 

dt dq i dq i dq i 


1 = 1 , 2 . 


(2.1) 


where KE and PE are kinetic and potential energies. q k ’s are the generalized coordinates and 
are defined by 

[q, q 2 1 = [0 i>] (2.2) 


GF’s are generalized forces. For the system without applied forces at the end effector, the 
generalized force vector is given as 

GF = [GF, GF 2 ] t = [T Of (2.3) 

where T is an applied torque at the joint. [ f represents the transpose of the matrix. 

The total kinetic energy has three parts, i.e., the kinetics energy of the arm (KE ,), 
the rotor of the actuator (KE r ), and the payload (££p. The mathematical expressions of these 
energies are 

KE. = - [ rJ r. dm (2.4) 

1 2 
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(2.5) 


KE. 


— Trace f r r r T . dm 

2 J rotor r r 


and 


KE=- Trace f 
p 7 J i 


rr‘ 

payload P P 


dm 


( 2 . 6 ) 


where f is the velocity of a differential mass, the subscripts Z , r, and p are for the link, the 
rotor, and the payload, and dm is the differential mass. Note that applying the Trace operator 
on the kinetic improves the computational efficiency for the rigid-body modeling since time- 
invariant terms can be separated from time-variant terms through the operator. 

The potential energy comes from strain energy and gravitational energy. The 
mathematical expressions are given as 


PE. 


-if 

7 J Hi 


El 

link a 


&V, 2 

(_J) dx 

ax 2 


(2.7) 


and 

PE = - f rj g dm - [ r T g dm (2.8) 

* J link 1 J payload P 


where EI a is the bending rigidity in the xy plane, g is a gravitational acceleration vector. 
According to the inertial coordinates in Figure 3, the gravitational acceleration is defined as 
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g = [0 0 -9.80] 7 (ml sec 2 ) 


(2.9) 


j££ and PE are expressed in terms of the generalized coordinates . The absolute positions and 

velocities in Equation (2.4)-(2.8) are formulated by the kinematics of the ERLS. [Ref. 1] 

Two sets of equation of motion for a single-link flexible manipulator are obtained 
from the derivation of above, which are nonlinear, coupled, second-order, ordinary differential 
equations represented as follows, 

m ee® + = fo (2J0) 

+ = /„ (2,U) 

where u represents the nodal displacement. m w and are effective masses for the large 

motion and small motion, and are coupled masses between large and small motion, g n 

is gyroscope for small motion, is stiffness for small motion, and / e and are 

load for large and small motion. The effective and coupled masses, the gyroscopic matrix, and 
load are nonlinear in 0 or 0 . In order to separate the applied torque from other terms in / e 

, let / 6 = h e T. [Ref. 3] In this study, a motion control is designed to control the tip position 

of the arm. With a small deflection assumption, the tip position can be approximated and 
represented by a total angle <p , where 9 = 0 + u/L. The control design therefore requires 
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a mathematical model for the total angle. Equations (2.10) and (2.11) can thus be rewritten into 
an explicit form for the applied torque as follows, 

^ee6 + w 8„ c = *0 + T (2,12) 

m n 0 ® + + W + *n u = 4 (2,13) 


Since 



Equations (2.12) and (2.13) can further be rearranged in terms of the total angle as 

w ee<P + ( m dr\ ~ = h 9 + T (2.14) 


m »8«P + K„ ‘ m J L )» + 8^ + *„ u = /, 


(2.15) 


From Equations (2.14) and (2.15), an uncoupled equation for the total angle <p is obtained by 
eliminating the small deflection u as follows, 

N <p + F c = T (2.16) 


where 
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N = m M - D 


*■.-*>(/, 6 «) “A. ( 217 > 

D = W *i ~ m ee/ L 

" m J L 

3. ELECTROHYDRAULIC ACTUATION 

The flexible arm is driven by an electrohydraulic actuator of which the dynamics 
is an integral part of the total system model. The dynamics of the electrohydraulic actuator 
include servovalve dynamics and hydraulic motor dynamics. A simplified description of 
servovalve dynamics was provided by MOOG, the manufacturer of the servovalve. A single 
equation which presents the dynamics is given as 

Q = KlJP v (2-18) 


where Q is the flow delivered from the servovalve, K is a valve sizing constant which 
contributes to hydraulic system damping. / is an input current, and P v is the valve pressure 

drop, i.e., P s - P L , where P g is the supply pressure and P L is the load pressure drop. 

Motor dynamics consists of a form of continuity equation and the torque output 
equation [Ref. 13]. They are written as follows 

Q = + V: + < 219 > 
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and 


T~H,r t D m (2.20) 

where D m is the motor displacement, C m P L is the leakage flow in the motor, (V f / 4p e )P L is 

the compressibility flow, and r| { is the torque efficiency. A detail account for the selection of 
hydraulic component for the system was included in [Ref. 3]. 

C. STATE SPACE REPRESENTATION 

Considering the plant dynamic model, Equations (2.16) through (2.20), these equations 
can further be rearranged as 

Mp + F e - x)P m P L (2-21) 

and 

Dj * C m P L * = K1 S [P W (2.22) 

Differentiating Equation (2.21) 

N<p + N# + F c = r\P m P L (2-23) 

Equation (2.21) and (2.22) can be rewritten as 
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(2.24) 





\Kl^Dj>-CJ> L ) 


\ * ) 



‘ 1 ' 


iAJ 


(Mp - F c ) 


(2.25) 


Substituting Equations (2.25) into Equation (2.24), and then substituting Equation (2.24) into 
Equation (2.23) give 

<p - v - u < 2 - 26 > 


where 



yp \ 


Mp + F c 


[WfJc /c . 4 MA. 2 q1 


K 4p 'S 

Vfi N Vp 


N K j 


Note that Equation (2.26) is a time-varying, nonlinear third-order ordinary differential equation 
represented in scalar form. 

In this study, the state space representation of Equation (2.26) will be used for 
the control action and the representation is given as follows 
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X = Bu +f 


( 2 . 27 ) 


where 


x = [ <p 4> <p ] r 

B = [ 0 0 btf 
u = I 

f = [ <P <P / 0 1 T 


The representation of Equation (2.27) will be used for the development in the control techniques. 
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ffl. CONTROLLER DESIGN (SLIDING MODE CONTROL) 


A. INTRODUCTION 

For many control system design problems, the designers do not have a detailed state- 
space model of the plant to be controlled, either because it is too complex, or because its 
dynamics are not completely understood. No nominal model should be considered complete 
without an assessment of its errors. These errors are refereed to as modeling uncertainties of the 
system. Therefore, a robust control design was attempted such that the control system will be 
insensitive to the modeling uncertainties. 

Based on our knowledge of physical mechanisms which cause differences between model 
and plant and our ability of representing these mechanisms, the representations of uncertainties 
may vary in terms of the structure they have. In practice, it is possible to represent these error 
in a highly structured parameterized form. These are usually the low frequency components. For 
a manipulator system, these parametric uncertainty or structured uncertainty may come from the 
imprecision on the manipulator mass properties, unknown loads, uncertainty on the load position 
in the end-effector, and inaccuracy on the torque constants of the actuators. However, there are 
always remaining higher frequency errors in the systems, which can not be covered in the 
parameter uncertainties. Also, these high frequency unmodeled dynamics is referred to as 
unstructured uncertainties. These unstructured uncertainties are usually caused by unmodeled 
structural modes, neglected time-delays in the actuators, or finite sampling rate. 

The greatest challenge of designing a robust controller is not only to minimize 
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performance sensitivity to uncertainties due to system parameters but also not to excite the high 
frequency unmodeled dynamics. The Sliding Mode Control has been known as a robust control 
that is able to work under the confined uncertainties of dynamic modeling error and 
environmental disturbances. The concepts of the SMC derived from the Variable Structure 
Control have been extensively studied in the Soviet Union for more than two decades. The SMC 
utilizes a high-speed switching control law to drive the plant’s state trajectory toward a specified 
surface (the sliding surface) and to maintain the plant’s state trajectory on this surface. 

The SMC using output models involves two filters, i.e., S dynamics and error dynamics. 
The 5 dynamics was designed to filter the uncertainties and the error dynamics was designed to 
obtain system error e. Figure 1 illustrated the sliding algorithm from the filter point of view, and 
also the control system block diagram is shown in Figure 4. Since a high-speed switched control 
about the sliding surface is not favorable to mechanical systems, a boundary layer thickness was 
introduced into the SMC such that a smoothed control is achieved. 

In this chapter, three sliding mode controls will be presented, which includes the sliding 
control with a first-order sliding condition, the straight sliding control, and the versatile sliding 
control. Before presenting the control algorithms, the uncertainties and the sliding surfaces will 
first be described. 
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Figure 4 Block Diagram of Sliding Mode Control 


SYSTEM EQUATIONS AND CONFINED UNCERTAINTIES 
From Chapter II, the dynamic equation of a flexible single-link arm was written as 

X = Bu + / ( 3J 

X = [Xj x 2 x 3 ] t 
u = [u] 

f - [/. / 2 /a ] r 


* - [ *1 *3 r 








Note that B and / are, in general, nonlinear functions of output variables and time. Equation 
(3.1) can be seen as a canonical form of physical models for nonlinear time-varying dynamical 
systems. 

To control system, a nominal mathematical model of the system can be obtained as, 

X = Bu + / ( 3 - 2 > 

B and f are the nominal values of B and /, which can be estimated from the theoretical 

evaluations or experiments. Due to parameter uncertainty, the discrepancies between the model 
and the physical plant are specified by A B and A f, where 

A B * B - B (3.3) 

A/ •/- / 

The sliding control assumes that uncertainties and disturbances are bounded. The 
uncertainties associated with the model are confined as 

|AB| s f, (3.4) 

|A/ I s; y 

where |*| denotes a norm of • which is a vector or a matrix. Note that in this study the norm 
of A , for instance, is defined as 

iai = v /Ps(a 5 x)l; < 3 - 5 > 
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where [eigC*)],^ stands for the maximum eigenvalue (i.e., the spectral radius) of • . This matrix 

norm is called the spectral norm and the corresponding vector norm is called Euclidean vector 
norm. The confining parameter p and y can be found as 



(3.6) 


where the components of ^ (or • min ) are the maximum (or minimum) values of the 

corresponding components of •. It is noted that p and y are non-negative numbers. It also 
implies that the system parameters and their nominal values are related as 


A 

B = 



^mn + ^aun 
2 

/ + f ■ 

J mix ■/ mm 


(3.7) 


The greatest value of confinement of A B is defined as 


(3.8) 
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C. SLIDING SURFACES 


In this single-link flexible arm study, sliding surface will only have one because of one 
control-input. The sliding surface is defined to be S = 0 where s is a sliding variable. The 
sliding variable directly relates the system error as 

S = R T e (3-9) 


where 

S = [ s t ] 

* * [ «1 e 2 ‘i f 

* - [ r i h r 3 f 

e is the tracking error, i.e., e = X - X d , where X d is a desired output vector, and is defined 
as 

= f X ld X 2d X 3d 3 

R is to locate the poles of the error dynamics on the sliding surface. S can be interpreted as an 
input to the error dynamics, and the error is a filtered version of S. 

An integral control can be introduced into sliding surface to eliminate steady-state error 
as follows, 
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S - * 7 «- C r /«* 


( 3 . 10 ) 


where 


G = [ «2 *3 f 


Note that once the poles for each sliding surface were assigned, the R and G are both constant 
vectors. 


D. SLIDING CONTROL WITH A FIRST-ORDER SLIDING CONDITION 

The relationship between control input u and S dynamics can then be obtained by 
differentiating Equation (3.10) 


S * R T i + G T e 


(3.11) 


Since 


X = Bu + / 
€ =X-X d 


Equation (3 11) can be rewritten as 

S = R\Bu * f - XJ+ G T e (3.12) 


The role of the control input u is to control the S dynamics such that the sliding surfaces can 
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be reached within a finite time. Once the sliding surface is reached, zero steady-state tracking 
can thus be achieved on the sliding surfaces provided the error dynamics is stable. 

The Lyapunov stability criterion is used to derive a sliding condition which specifies the 

S dynamics. A Lyapunov function V is picked such that V = — S T S. By applying the Lyapunov 

2 

stability criteria, ^o- a sliding condition is written as 

S T S z 0 (3- 13 ) 

The sliding condition assures the attractiveness of the S dynamics toward the sliding surfaces. 
To be able to adjust the sliding speed 5 , Equation (3.13) is rewritten as 

S r Ss-q|S| ( 314 > 

where ^ is a non-negative real number. It is worth to note that Equation (3.14) is a first-order 
sliding condition. 

The sliding condition specifies a desired dynamics of S on which S slides toward the 
sliding surfaces and the steady-state error can be eliminated. To obtain the desired S dynamics, 
the control input can be designed using predictor-corrector scheme as 

u = (R t BY\& - k sgn(S)) ( 315 ) 

where £, i.e., nominal control input, can be obtained by letting 5 = 0 - which will give 
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& = R T (X d -h - G T e 


(3.16) 


Also, 


sgn(S)= [ ] 


1 s^O 

-1 s ,<0 


Note that (R t b) is a scalar- Also, note that (R T By l & is a predictor and ( R T B) l k sgn(S) is a 

corrector. The gain it is a non-negative real number and is determined by using matrix-norm 
techniques. By substituting Equation (3.12) into Equation (3.14), an inequality will be given as 

S T [R\Bu+f-X d )+G T e] s -t||S| (3-17) 

The control input (Equation (3.15)) can then be substituted into Equation (3.17), which will give 
S*[(Jl T B)(R t B)'\H- k sgn(S))+R T (f-X d )+G T ey-T]lS\\ (3.18) 

Also, Equation (3.16) will give 

R T X d = a + R T f - G T e 

Substituting into Equation (3.18) and rearranging give 

S T {R T (f-f) +[(R T B)(R T BY l ~l]u - (R t B)(R t B) l )k sgn(S)}^-q |S| ( 319 ) 

The uncertainties described in Equation (3.3) are applied to Equation (3.19) and the sliding 
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condition is given as 


S T \R T Af+(R T AB)(R T By'&-[l+(R T AB)(R T By 1 ^: sgn(S)}*-r\ ||S| ( 3 - 2 °) 


In order to quantity gain k, the matrix-norm technique is applied to Equation (3.20) and 


S t {R T Af+(R t AB)(R t B) 1 H -[/+(* t AB)(R r i)- , ]jfc sgn(S)} 
z\S[{\R\\(\Af\ + \AB[\(R T B) l fi[)-k(.l-\RABmR T By l sgn(m} 


(3.21) 


Note that 


-kS T (R T AB)(R T B) l sgn(S) s k\S\\RAB\\(R T B) 'sgn(S)\ 


~kS T sgn(S) s -*»S| 


and 




Therefore, 


\S\{\RK\Af\+\ABlM T B) l fi\)-k(l-\RAB\M T By'sgn(S)\}± T)\S\ < 3 - 22 ) 


Because of the uncertainty confinements (Equations (3.4) and (3.6)), Equation (3.22) becomes 
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(3.23) 


Finally, the gain k is found as 

^ |Jg»(Y + Pi(* r l)- l fi|)+Ti (3 24) 

1 - |J?A B 1l(lt T B)' 1 sgn (S) || 

As stated before, the gain £ is a non-negative number. Thus, a sufficient condition is given 

1 > |JtA£||(Jt r i)- 1 S£it(S)| ( 3 - 25 ) 

The Lyapunov stability is guaranteed as long as the gain jfc is chosen according to Equation 
(3.24). For control design purposes, the minimum value of k is selected since the least control 
effort is desired. Because of the discontinuity of the sgrt functions, the control is called switched 
control and causes chattering, which is unfavorable to the mechanical devices. Equation (3.15) 
can be substituted into Equation (3.12) to obtain the S dynamics as follows 

S+k(R T B)(R T By l sgn(S ) =R T (f-f)-(R r A B)(R r B) 1 u ( 3 - 26 ) 

It is observed that the right-hand-side of Equation (3.26) consists of uncertainties, error, and 
desired trajectory dynamics, which are the excitations to the S dynamics and are to be filtered 
out by the first-order S dynamics. 

To smooth the control law, a boundary layer with thickness $ is introduced where the 
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thickness is a real positive number. Outside the boundary, |s, | ;>4>, the control law u is designed 

to satisfy the sliding condition, Equation (3.14), which guarantees boundary layer attractiveness. 
Inside the boundary, |s,|«{>, the control law will impose a smoothing process to the S 

dynamics. The general law is written as 

u = - k (3-27) 

where 

H = R T (X d - f) - G T e 


The sat function is defined as 

sat(^-) = sat (= • 

9 9 s t 

♦ 


ki l<4> 


The S dynamics of Equation (3.26) is smoothed within the boundary layer, i.e., 

S+-(R t B)(R t B) 1 S=R T (f-f)-(R r A B)(R T B) x u (3.28) 

4 > 


Equation (3.28) represents a first-order low-pass filter. The boundary layer thickness <j> 
determines the response speed of the S dynamics and the values of <{> can be either constant or 
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time-varying. 


To obtain a time-varying thickness, a maximum bandwidth A. m , T is assigned to the S 

dynamics within the boundary layer (Equation (3.28)). The analytic properties of matrix-norm 
are to be used to quantify the bandwidth. The spectral radius of the bandwidth is defined as 

X = -(R t B)(R t B ) 1 

4> 


Thus, the boundary layer thickness is designed such that Also, 


Therefore, 


~(R T B)(R T ByU - ^R T B)(R T B) l \ 
<t> <t> 


(3.29) 


Is - \iR T B)(R T B) l \ 
<l> 


(3.30) 


Since 


M T B)(R T By^ i l + \\R£iBn(R T B)- l \ ( 3 - 31 ) 


and by the confinements of uncertainty (Equation (3.3)), the spectral radius is then found as 

xAl + mBUiRTBy'D 

9 
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Hence, the maximum bandwidth can be assigned as 


<p 

Therefore, the thickness <}> is 

4> - ~^-(U\\R£iB\\(R T By l \) (3.32) 

In the Equation (3.32), the maximum bandwidth determines the boundary layer thickness^ 

which controls the tracking accuracy and response speed. By selecting a proper value of , 

Equation (3.32) can assure that bandwidth will never exceed . 

In summary, the sliding control with first-order sliding condition is designed to have both 
S dynamics and error dynamics low-pass filters where an integral error can be added on error 
dynamics to eliminate the steady-state error. The S dynamics is a first-order low-pass filter, in 
which the uncertainties are filtered. A detailed synthesis view of S dynamics of sliding control 
with a first-order sliding condition is shown in Figure 5. The continuous control of Equation 
(3.26) can be used to replace the switching control and eliminates the chattering. A detailed 
block diagram of sliding control system is illustrated in Figure 6. 
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Figure 5 A Detailed Synthesis View of Sliding Mode Control with 
First-Order Sliding Condition 


R X«j 



Figure 6 A Detailed Block Diagram of Sliding Mode Control with 
First-Order Sliding Condition 
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E. STRAIGHT SLIDING CONTROL WITH FIRST-ORDER PLUS INTEGRAL 
SLIDING CONDITION 

The straight sliding control begins with the same sliding surface which was defined in 
Equation (3.9). By differentiating Equation (3.9), the relationship between control input and S 
dynamics can be found 

S = R\Bu * / - X} (3.33) 


The Lyapunov stability criterion will be used to derive straight sliding condition and 
specify the desired S dynamics such that the sliding condition can be reached in a finite time. 
A Lyapunov function is defined as 




(3.34) 


By applying the Lyapunov stability criterion, y s q» the sliding condition is obtained 


S T (S+u*f'sdt) a 0 


(3.35) 


Equation (3.35) indicates that a first-order plus integral sliding condition representing a stability 

criterion for an equivalent mass-spring-damper system in which dt is a equivalent 

Jo 

displacement. The g> w 2 JJs dt provides an additional restoring effort and the dynamic behavior 
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can thus be tuned by the equivalent spring constant . Furthermore, Equation (3.35) can be 
written as 

S r (S+u 2 f's dr)s-Tt|S|| (3.36) 

" JO 

where n is a positive real number. The sliding condition, Equation (3.36), describes a desired 
S dynamics on which S slides toward the sliding surface. 

In order to obtain the desired S dynamics, the control input can be obtained as 

u = (R T B)-\it - ksgn(S)) ( 3 - 37 > 

where £ is obtained by letting S+u> 2 jSdt=Q, 

a = R t Bu = R T (X d -f)-u n 2 jsdt (3.38) 

The gain k is determined by substituting Equation (3.33) into Equation (3.36), which gives 

S ^R t (Bu +f-X d )+u> 2 fSdt)s -T) |S| (3.39) 

Equation (3.39) can be further simplified by using Equation (3.37) and Equation (3.38), 

S T {R T bf+(R T £iB)(R T By'u-[l+(R T bB)(R T By^ sgn(S)}z - n ||S| ( 3 - 40 ) 
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It is noted that Equation (3.40) has exactly the same form as Equation (3.20) has. Thus, the 
same matrix-norm techniques used in Equations (3.20-25) will be utilized to quantify the gain 
k, which gives 


hi 


l*|(y+p|(Jt , 3)- t *|) 

l-lRASfi!(tf r i)-Vi(S)! 


(3.41) 


The S dynamics can be formed in terms of uncertainties by substituting Equations (3.37- 
38) into Equation (3.33), which gives 


S+k(R T B)(R T B)~ l sgn(S) +(R T B)(R T Br l <j> m 2 fsdt 
=[(R T B)(R T By l -l]R T X d +R T f-(R T B)(R T B)~ x R T ? 


The right-hand-side of Equation (3.42) represents the excitations to the S dynamics, which 
consist of the desired trajectory dynamics, uncertainties, and disturbances. In other words, the 
5 dynamics can be treated as a filtered version of excitations and filter out most unwanted 
uncertainties and disturbances. The sgn function causes chattering which will produce 
undesirable noise and mechanical wear. A continuous control is then developed to eliminated 
chattering. 

A boundary layer with thickness <J> is introduced to smooth out S dynamics. The control 
algorithm is defined as 

« = (J? r B)-‘^ - k «tf(|)j 

where 
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a = R\x d -f)- u .'rs it 


sct(^) <= saK-i) 
<P 4> . 


s. 

<t> 


sj*4> 


KM> 


Thus, a smoothed first-order low-pass S dynamics can be obtained as 


S+-^(R T B)(R T By l S+(R t B)(R t B) - 1 u 2 J Sdt 
{{R t B)(R t B) 1 T X d +R 7- (R T B) (R t B)'R T f 


(3.43) 


The S integral introduces an integral control to the S dynamics and guarantees a zero steady- 
state values of S which will drive the steady-state error dynamics to zero. 

From Equation (3.43), the maximum bandwidth A m . T is assigned to S dynamics within 

the boundary layer. The analytic properties of matrix-norm are to be used to quantify the 
bandwidth. The bandwidth (A.) is defined as 


k 2 = (R T B)(R T By x u 2 


Therefore, 


(r t b)(r t b)-'u 2 s Yr t bxr t b )-% 2 


(3.44) 
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X 2 s \tR T B)(R T By l \u m 2 


(3.45) 


Since, 


|(/J T B)(U T i)* 1 S S l + lR*B\M T B) l [ ( 3 46 > 


The bandwidth is then found as 

X 2 s (l + U?AJMi(il r ir 1 |)(o. 2 (3-47) 


Hence, the maximum bandwidth is 

X^ 2 = (1+{J?A£| |(U t J5) _1 Dm,, 2 < 3 - 48 ) 

Since the maximum bandwidth has been chosen such that the bandwidth of the S dynamics will 
never exceeds X m „ and the unmodeled dynamics can be filtered. From the damping of Equation 

(3.43) and also J? = AJJ + B 

~(R t B)(R t £)-* = ~\(R T &B)(R T B)+l} 


Also, 


—[(R T A.B)(R T By l ♦ l] i -(R T &B)(R T B) 1 

4 > 4 > 


Therefore, a lower bound of damping is 


(3.49) 
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(3.50) 


2£,X, = ±(R t AB)(R t B) 1 
4 > 

Note that (R t &B) and (R T B) are both scalar. Thus, 

(tf r A£)(rt r Br 1 

X, is replaced by , which gives 






\(RAB)(R t B)' 1 


Hence, a more conservative lower bound of damping is obtained as 




2^4*J 


\(RAB)(R r JB)" 1 


Therefore, the thickness $ is. 


4 > = 


( k \ 




(RAB)(R T B) 


Tfi \-1 


(3.51) 


(3.52) 


(3.53) 


X m>T is given by Equation (3.48), where the spring constant is the only control parameter 
for the bandwidth. Once is determined, another tuning parameter ( L can be properly 
adjusted. It is shown in Equation (3.53) that also governs the activeness of S dynamics where 
the activeness refers to the magnitude level of S response. The greater C L is, the thinner the 
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boundary layer thickness $ is, and the less active of S dynamics is. Also, the less active S 
dynamics has less influence over to the error dynamics. However, the thinner boundary layer 
may require higher sampling ratio to eliminate the chattering. 

In summary, the straight sliding control is designed to have two lower-pass filters in 
which S dynamics is designed as a first-order low-pass filter with an integral of S . Figure 7 
shows a detailed synthesis view of straight sliding control system. The integral of S provides an 
ability to make the steady-state value of 5 to zero, and the steady-state error will be driven to 
zero. The straight sliding control provides two independent tuning parameters to adjust the S 
dynamics. The tuning parameter of spring constant (u ( 2 ) provides a restriction to the bandwidth 

of S dynamics which the unstructured uncertainties (unmodeled dynamics) will be rejected. In 
addition, the C L gives a fine tune on the thickness of boundary layer. Figure 8 show a detailed 

block diagram of straight sliding control system. 

F. VERSATILE SLIDING CONTROL WITH SECOND-ORDER SLIDING 
CONDITION 

The idea of versatile sliding control is to use a low-pass S dynamics filter to filter out 
unwanted high-frequency noises that include unmodeled dynamics and uncertainties, and then 
followed by a high-pass (or band-pass) error dynamics filter. 

In order to make error dynamics a high-pass filter, an additional zero is placed into the 
error dynamics. Thus, the error dynamics is modified to be 
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Figure 8 A Detailed Block Diagram of Straight Sliding Control 
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S +zS =R T e+G T J"edt 


(3.54) 


z is a scalar where a zero is assigned for sliding surface. The relationship between input control 
and S dynamics can be found by differentiating Equation (3.54), which gives 

S+zS =R T (Bu+f-XJ)+G T e (3.55) 


Now, Equation (3.55) is a second-order form of sliding condition, which specifies S dynamics 
in order to reach sliding surface. Thus a Lyapunov function to suit the second-order sliding 
condition is 

V = —S T S+—S r w 2 S (3.56) 

2 2 


By differentiating Equation (3.56) and imposing sliding speed parameter, the Lyapunov stability 
criteria becomes 

s\s+v>*S)z-t\ |S| (3-57) 

where sliding speed, t], is a positive real number. 

According to the second-order sliding condition, the switched control law will be 

u = ( R t B)-'( 6 - k sgn(S )) (3-58) 

where the nominal a can be obtained by letting s + u n 2 S = 0 
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£ =R T (X,-f)-G T e+ z S-u> 2 S 


( 3 . 59 ) 


and the gain k can be found based on the same knowledge as before. Substituting Equation 
(3.55) into Equation (3.57) and further simplifying by Equation (3.58) and Equation (3.59) give, 

S T {* T Af+(R T AB)(R T i)-‘<2-[l +(J? t AB)(R T J?r']* sgn(S)}<;-n JS|| ( 3 - 60 ) 


Applying the matrix-norm technique, the gain k can be found as 


\R\\{y+V[(R T B)- l u\)+r\ 


l-\RAB\\\\(R T By l sgn(S)'! i 


(3.61) 


Because of chattering, the switched control law is modified to be a continuous control law as, 


u = (tf r 2f)- 1 |fi - * sol(^) 


(3.62) 


Outside the boundary layer, the control law is governed by Equation (3.58). While inside the 
boundary layer, Equation (3.62) governs the control law. Thus, the definition of sat function is 
defined as 


sati^-) = sat(-±) 

<P <P 


i, 

'gn(-j-) 

9 




k, l<4> 
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The continuous S dynamics can then be obtained by substituting Equation (3.59) and Equation 
(3.62) into Equation (3.55) 

S+\±(R t B)(R t B)-' +{l-(R T B)(R T By 1 ]As+(R T B)(R T B)-' u h 2 s 

<J> J (3.63) 

=R T f-(R T B)(R T B)- 1 R T f +[(R T B)(R T By l -l]R r X d +[l-(R t B)(R T By l }G T e 


To find time-varying boundary layer thickness, the maximum bandwidth is first to be evaluated 
as in the section D (straight sliding control), 


* * = (1+ ^RAB\ \(R t B)' 1 Du. 2 


(3.64) 


The damping of S dynamics is 


4(R 7 "fi)(K r B)' 1 + (f-(. R T B)(R T Sy l )z 
<»> 

= ~{(R T AB)(R t B)" 1 +l) - (R T AB)(R T By'z 
<P 


(3.65) 


Also, 


~{(R T AB)(R T By 1 +l) - (R T AB)(R T B) 'z 
<P 

* ^-(R T AB)(R T B)-' - (R T AB)(R T B) l z 
4 > 


(3.66) 


Therefore, a lower bound of damping is obtained 
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(3.67) 


2C,X| = -(R T AB)(R T By l - ( R T ABXR T B) l z 
<t> 


and 



(R T ±B)(R T By l 

2k, 



(3.68) 


X, is replaced by X m . T , which gives 


C/ 


' (R t HB)(R t B)- 1 



Thus, a more conservative lower bound of damping is obtained as 



(* r Al?)(* r fi)- 1 



Therefore, a time-varying boundary layer thickness is, 


kjR t &B)(R t B)~ 1 

2Z L k nix HR T *B){R T Br l z 


(3.69) 


(3.70) 


The relationship of tuning parameters (C z , ^ , z) are evaluated in Equation (3.70). The 

bandwidth for the S dynamics is given by Equation (3.64), in which the spring constant is 

the only control parameter for the bandwidth. Equation (3.70) provides two other control 
parameters (C L , z) which give more capability of tuning and shaping the S dynamics. 
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In summary, the versatile sliding control provides a second-order low-pass filter (S 
dynamics) and a high-pass (or band-pass) filter (error dynamics). Figure 9 illustrates a detailed 
synthesis view of versatile sliding control system. The 5 dynamics as low-pass filter is designed 
to reject the high-frequency uncertainties (unmodeled dynamics, error dynamics and 
disturbances). Low-frequency uncertainties are then filtered by high-pass (or band-pass) filter. 
A detailed block diagram of versatile sliding control system is shown in Figure 10. 
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Figure 9 A Detailed Synthesis View of Versatile Sliding Control 



Figure 10 A Detailed Block Diagram of Versatile Sliding Control 








IV. RESULTS 


To perform simulation and implementation, the structured uncertainties were first to be 
determined. The nominal values of the system parameters and their bounds were estimated. A 
non-robust controller was designed to estimate the nominal values and the bounds [Ref. 14]. The 
time history of these two parameters for two different payload (0 kg and 0.85 kg) illustrated in 
Figure 11 (a)- 11(g). Thus, the value for each parameter can be determined and £, p , f, andy 

can be found in the Table 2. 


Table 2 


Payload 

P 

A 

B 

Y 

/ 

0.00 kg 

295.00 

[0 0 11005.00f 

677177.18 

[0.9022 10.7036 -57732.00] T 

0.85 kg 

218.95 

[0 0 5932.20f 

35684.40 

[0.8557 5.6229 -35097.00] T 


For conservative reasons, set of parameters for the 0.00 kg case was selected and used 
in the simulation and implementation throughout the whole research. The computer simulation 
of the control system was performed on a 16 Mhz 80386 IBM compatible personal computer 
with MATLAB programs. The MATLAB simulation programs and the FORTRAN 
implementation programs are listed in the appendix A and B. 

The SMC control algorithms used for simulation and implementation on the single-link 
flexible arm are the sliding control with first-order sliding condition, The straight sliding control, 
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and the versatile sliding control. The following sections will evaluate the performances of these 
three sliding controllers. To show the robustness to the uncertainties due to varying payload and 
to perform the tuning procedure such that the unmodeled high-frequency will be rejected are the 
main goal of this research. The control system performance will be evaluated according to the 
tracking accuracy, the response speed, and the overshoot. 

A. RESULTS ON SLIDING CONTROL WITH FIRST-ORDER SLIDING CONDITION 
PLUS INTEGRAL ERROR 

The system was simulated with a time interval of 0.001 (sec). The controlling parameter^ 
for sliding speed was selected to be 10. An unit-step response of continuous sliding control 
system of simulation and experiment was first investigated. Figures 12(a)-12(d) show the system 
responses (i.e., controlled tip position (<p), control input, S response, and boundary layer 
thickness (<j>)) of simulation without payload where = 10 (rad/sec). The poles of error 

dynamics were at p=[-10 -10 -10]. The tip position response shows large overshoot. Time- 
varying boundary layer thickness ($) was computed to achieve the continuous sliding control. 

It was observed that |s, | < <J> and the S dynamics was governed by a first-order dynamics where 

the controlling parameter X m , T determines the response speed. The large overshoot will be 

minimized by fine tuning the S dynamics such that the system response speed was primarily 
determined by placing the poles of error dynamics. The tightness of the boundary layer can also 
improve system response speed. 

A higher bandwidth was selected, i.e., A. mjT = 500 (rad/sec) to reduce the boundary 
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layer thickness and to reduce the overshoot. Figures 13(a)-13(d) present the performance with 
no steady-state error and the overshoot was improved. It should be noted that the greater 
bandwidth will allow high-frequency unmodeled dynamics to excite S dynamics. It was observed 
that with a greater k mMT , the response speed is faster. However, the level of control input is 

increased. With the tighter boundary layer, the level of overshoot was improved. Because 
boundary layer thickness becomes thinner and the sampling rate was not high enough, the 
occurrence of chattering is not surprising and the current (i.e., control input) saturation occurs. 
By selecting the higher sampling rate to eliminate chattering and saturation, i.e., A t = 0.0005 
(sec), Figures 14(a)-14(d) show the further improvement of control performance in the 
simulation without payload. Figures 14(e)-14(h) show the experimental results. Figures 15(a)- 
15(d) show the responses in simulation and Figures 15(e)-15(h) show those in experiment with 
payload (0.85 kg) when A. m , T = 500 (rad/sec) and simulation time interval At = 0.0005 (sec). 

The next test was to examine the tracking performance. The desired trajectory was 
defined as 



t 

1 

3-r 

0.6 


for 0 <; f < 1 
for 1 s, t < 2 
for 2 s t < 2.4 
for t > 2.4 


(4-1) 


Figures 16(a)-16(h) show the tracking performance of the system in simulation and experiment 


without payload and Figures 17(a)-17(h) show the tracking performance with payload (0.85 kg). 


A greater bandwidth (A. m>T = 1700 rad/sec) was selected to tight the boundary layer. 


With the trajectory control, the overshoot can be further minimized. The improvement 
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on reducing the overshoots is very important for robot applications. The overshoot in the 
experimental results reflected the simulation prediction. The time history of the boundary layer 
thickness in the simulation differs from that in the experiment. This difference can be neglected 
since the change of the boundary layer thickness compared to its magnitude was quite small. The 
controller is robust to varying payload because the system responses were not affected by 
changing the payload. 

In summary, the maximum bandwidth is the only parameter to reject the unmodeled high- 
frequency dynamics and at the same time to control tracking accuracy in this control algorithm. 
Therefore, a trade-off exists between the robustness to the unmodeled dynamics and the tracking 
accuracy. With an integral control in error dynamics, the steady-state error was eliminated. The 
controller is indeed robust to the uncertainties due to varying payload. 

B. RESULTS ON STRAIGHT SLIDING CONTROL 

Step response will first be studied to examine the performance of the straight sliding 
control. An equivalent spring constant and a damping ratio were selected as = 10 (rad/sec) 

and C L = 1 simulated at a t = 0.0002 (sec), while the poles of error dynamics were selected 

as p=[-10 -10]. Figures 18(a)-18(d) show results of a step response, the control input, and the 
S and of the simulation without payload. The control presents not only overshoot but also 

long settling time. Thus, damping ratio (C L ) will be used to fine tune S dynamics such that the 
overshoot will be minimized and the settling time will be improved. The equivalent spring 
constant (co b 2 ) will be lowered to reject high-frequency excitation of the S dynamics. 
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The effects due to the damping ratio { L will be examined. The damping ratio = 10 
was selected to reduce the boundary layer thickness while a* = 10 (rad/sec). Figures 19(a)- 
19(d) show the performance of step response, the control input, and the S and $ of the 
simulation without payload. Figures 19(e)-19(h) show the step response of experiment at C* = 

10 and = 10 (rad/sec). The overshoot was minimized and the response speed was faster than 
before. 

The next test is to lower the i.e., lowering the bandwidth of the S dynamics to 

reject the high-frequency excitations, a* = 1 (rad/sec) and ( L = 10 were selected. Figures 

20(a)-20(d) present the control performance (without payload) where the step response has 
steady-state error. To reject the unwanted high-frequency excitations, the control system has to 
pay the price. However, the damping ratio C L can be increased to remedy the sluggish of the 

response. Figures 2 l(a)-21(h) show the control performance of = 1 (rad/sec) and C L = 100 

for both simulation and experiment without payload. Figures 22(a)-22(h) present the control 
performance of = 1 (rad/sec) and = 100 for both simulation and experiment with 

payload (0.85 kg). It is noted that the control performance showed no difference in changing the 
payload. 

The last investigation was to examine robustness to varying payload of the tracking 
performance. The desired trajectory was designed as in Equation (4.1). Figures 23(a)-23(d) show 
the tracking performance for the simulation at = 1 (rad/sec) and C t = 100 without payload. 
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Figures 23(e)-23(h) show the tracking performance for the experiment at = 1 (rad/sec) and ( L 

= 100 without payload. Figures 24(a)-24(h) illustrate the tracking performance for the 
simulation and experiment at o || = 1 (rad/sec) and £ L = 100 with payload (0.85 kg). Theta,, 

= 1 (rad/sec) was selected not to excite the unmodeled high-frequency dynamics. 

Note that the control performance was not affected by the payload. The greater the 
damping ratio C L , the better the tracking performance and the greater the control effort is 

needed. Although the smoothed control law was used to compute the current required to the 
actuator, a little chattering occurred in the experimental results while the arm was trying to stand 
still. 

In summary, the straight sliding control utilizes the first-order plus integral sliding 
condition to provide tracking accuracy, disturbance rejection, and stability. The S dynamics 
plays a key role to obtain a desired system response since the error dynamics is affected by the 
output of S dynamics. The two control parameters £ z and o^ 2 are adjusted to accomplish the 

control task without conflicting. The damping ratio £ L was selected to tight the S dynamics 

performances, while the q r 2 in the sliding condition to keep the excitation from entering into 

the S dynamics. Therefore, the trade-off between tracking accuracy and robustness to 
uncertainties (due to unmodeled high-frequency dynamics) has no longer exist. Also, the 
robustness to the varying payload was proved by the simulation and the experiment. 
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C. RESULTS ON VERSATILE SLIDING CONTROL 

The step response of versatile sliding control was to test the control performances. To 

keep the bandwidth of S dynamics low, the o ^ 2 was selected to be 1. The zero ( z ) not only 

provides the filtering level of error dynamics but also affects the boundary layer thickness. Thus, 
the z should be chosen carefully such that the S dynamics can be fine tuned. For a given value 

of C £ , the smaller z , the larger boundary layer thickness. Also, greater boundary layer 

thickness gives more active S dynamics. Given z = 1 and f £ = 10, Figures 25(a)-25(d) show 

the step responses, the control input, and the 5 and <j> of simulation, while the poles of error 

dynamics were p=[-10 -10 -10]. Figures 25(e)-25(h) show the experiment results. The overshoot 
exhibited in the control motion. With the same C £ » Figures 26(a)-26(d) present the responses 

of a greater value z = 10 without payload. The overshoot was reduced by reducing the 
boundary layer thickness. However, the steady-state error occurred. Figures 26(e)-26(h) present 
the results for the experiment without payload at z - 10 and f £ = 10. Figures 27(a)-27(h) 

present the simulation and the experiment with payload (0.85 kg) at z = 10 and C £ = 10. 

The next test was to examine tracking performance. The trajectory was designed as in 
Equation (4.1). Given z = 10, C £ = 10, and u ^ 2 = 1 (rad/sec), Figures 28(a)-28(h) show the 

performances of simulation and experiment. 

The overall tendency of the control system in the experiment has the same results as 
expected in the simulation. The smaller the boundary layer, the smaller the steady-state error, 
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and the greater control input is needed. The little chattering occurred in the straight sliding 
control was seen in the performance of versatile sliding control. It is noted that the robust to the 
varying payload is achieved by increasing the payload while the motion control will not be 
changed. 

In summary, the versatile sliding control provides three parameters to accomplish the 
control motion. The damping ratio and z were chosen to tight the boundary layer thickness 

and also the was set low to reject the excitations entering the error dynamics. The z not 

only gives the tuning on the boundary layer thickness but also make the error dynamics a band¬ 
pass (or high-pass) filter. Actually, the versatile sliding control provides more tuning parameters 
which make the trade-off between tracking accuracy and robustness to uncertainties (unmodeled 
high-frequency dynamics) no longer exist. Also, the versatile sliding controller is proved to be 
robust to the uncertainties due to varying payload. 
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Figure ll(a)-(d) The Time History of System Parameters 
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Figure 12(a)-(d) The Simulation Performance of SMC with First-Order Sliding Condition 
= 10 rad/sec, Ar = 0.001 (see), payload = 0.00 kg) 
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Figure 13(a)-(d) The Simulation Performance of SMC with First-Order Sliding Condition 
(L = 500 rad/sec, Ar = 0.001 (sec), payload = 0.00 kg) 
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Figure 15(a) 
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Figure 15(a)-(d) The Simulation Performance of SMC with First-Order Sliding Condition 
m , T = 500 rad/sec, At = 0.0005 sec, payload = 0.85 kg) 
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Hgurc !6(a)-(d) The Simulation of Tracking Performance of SMC with First-Order Sliding 
Condition (X^ = 1700 rad/sec, A t = 0.0005 (see), payload = 0.00 kg) 
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Figure 16(e)-(h) The Experiment of tracking Performance of SMC with First-Order Sliding 
Condition (A. = 500 rad/sec, payload = 0.00 kg) 
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figure 17(e)-(li) Tlie experiment of tracking Performance of SMC with First-Order Sliding 
Condition (A m „ = 1700 rad/scc, payload = 0.00 kg) 


63 






Figure 18(a)_ in _ Figure 18(b' 



(puy) S 


Figure 18(a)-(d) The Simulation Performance of Straight Sliding Control 
(U" = 10 rad/sec, = 1, Ar = 0.0002 (sec), payload = 0.00 kg) 
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Figure 19(a)-(d) The Simulation Performance of Straight Sliding Control 
(u # = 10 rad/sec, = 10, Af = 0.0002 (sec), payload = 0.00 kg) 
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Figure 19(e)-(h) The Experiment Performance of Straight Sliding Control 
(w, ~ 10 rad/sec, ( L = 10, payload = 0.00 kg) 
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Figure 21 (a)-(d) The Simulation Performance of Straight Sliding Control 
( = 1 ra d /s ec, = 100, Ar = 0.0002 (sec), payload = 0.00 kg) 


68 


Time (Sec) Time (Sec) 
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Figure 22(a)-(d) The Simulation Performance of Straight Sliding Control 
(<*>„ = 1 rad/sec, = 100, Af - 0.0002 (sec), payload = 0.85 kg) 
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Figure 24(a)_ _ Figure 24(b) 
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Figure 24(a)-(d) The Tracking Performance in Simulation of Straight Sliding Control 
(<*)„ = ! rad/sec, C L = 100, At = 0.0002 (sec), payload = 0.85 kg) 
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Figure 26(a)-(d) The Simulation Performance of Versatile Sliding Control 
(«„ = 1 rad/sec > z = 10, C £ = 10, Ar = 0.0002 (sec), payload - 0.00 kg) 
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Figure 27(a)-(d) The Simulation Performance of Versatile Sliding Control 
(<*>„ = 1 rad/sec, z = 10, C £ = 10, At = 0.0002 (sec), payload = 0.85 kg) 
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V. CONCLUSIONS AND RECOMMENDATIONS 


A. CONCLUSIONS 

The purpose of doing this research was motivated by the need of using a flexible 
manipulator in the modem robotic application and of designing a simple robust control (sliding 
mode control) algorithm enabling a single-link flexible manipulator to perform accurate tracking 
under the environmental disturbances and uncertainties due to varying payload. 

Because of the simple control structure of sliding mode control, the robust control system 
was implemented on a low-cost IBM-AT micro-computer for the flexible arm system. Unlike 
the high-speed switched sliding control, the predictor-corrector continuous control law was 
achieved by introducing a time-varying boundary layer. Three different forms of Lyapunov 
stability criterion (or sliding condition) were utilized not only to guarantee the stability of the 
control system but to provide different tuning capability. 

In the sliding control with first-order sliding condition, the only control parameter 

not only provides the tightness of the boundary layer but also gives bandwidth of the S 
dynamics. The tighter the boundary layer, the faster the response speed and the smaller the 
steady-state error. However, the trade-off between tracking accuracy and the robustness to the 
unmodeled high-frequency dynamics exists. The integral error control will drive the steady-state 
error to zero. 

The straight sliding control with first-order plus integral sliding condition provides control 
tuning parameters (<»>„ and £ L ) such that the rejection of uncertainties and tracking accuracy can 
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be achieved at the same time. The integral of S provides an ability to have the steady-state 
values of S zero for constant input to the S dynamics, and the steady-state error of the system 
will be driven to zero. The equivalent spring constant (q„ 2 ) provides an rejection to the 

unmodeled high-frequency dynamics. The lower bound of damping ratio (C L ) provides a fine 

tune on the thickness of boundary layer. The thinner the boundary layer thickness, the faster the 
response speed and the smaller the steady-state error. However, from the experiment, even with 
continuous control, the small chattering occurred while the arm wanted to stand still. This 
phenomenon has not been clear. 

The versatile sliding control algorithm provides a second-order sliding condition. The 
bandwidth of 5 dynamics was selected to be low so that the high-frequency unmodeled 
uncertainties will be rejected and the error dynamics will not be excited. The lower bound of 
damping and z in the 5 dynamics can then be fine tuning to tight the boundary layer so that the 
tracking accuracy can then be achieved. Therefore, the trade-off between the tracking accuracy 
and the robustness to the uncertainties has no longer existed. However, the same chattering in 
control input as in the straight sliding control occurred when the control was implemented. 

In summary, despite the use of simplified model, the simple structure sliding control is 
indeed robustness to uncertainties due to varying payload. Also, by providing more tuning 
capability, the trade-off between tracking accuracy and the robust to the high-frequency 
dynamics was released. 
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B. RECOMMENDATIONS 

The further works are recommended as follows : 

(1) To develop a systematic procedure to tune controlling parameters. 

(2) To clarify the small chattering occurred in the implementation while using straight 
and versatile sliding control. 

(3) To develop the multiple-link control system using ERLS. 
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APPENDIX A 


MATLAB SIMULATION PROGRAM 


DEFINATTONS OF THE PARAMETERS : 

* lambda = the max. bandwidth h = time interval 

* w = tunning parmeter to tunne the max. bandwidth 

* zeta = tunning parameter to provide more tunning capability 


= final time ml = payload 

= pole-placement of poles speed = sliding speed 

= B hat(uncrtainties) DB = DELTA B 

= fhat(uncertainties) Df = DELTA f 

= length of the arm rho = density of arm 

= cross asection area of arm E = modulusof elasticity 

= Area moment of inertia of flexible arm 
= deformation matrix differentiation w.r.t. deflection 


* 

g 

* 

K 

* 

Dm 

* 

Vt 

* 

Ps 

* 

Ctm 

* 

eff 

* 

betae 

+ 

Ip 

* 

Ir 

* 

betal 

* 

beta2 

♦ 

Cl 

* 

C2 

* 

nlm 

* 

nlnlm 

* 

lxnlm 

* 

swwsm 

* 

klei 

* 

maa 

* 

mab 

* 

mba 

* 

mbb 

* 

kb 

* 

fa 

* 

fb 

* 

theta 


* 


= gravitational acceleration vector 
= servovalve sizing constant 
= actuactor displacement 

= total compress volumeincluding actuactor lines and chambers 
= hydraulic supply pressure 
= total leaking coefficient 
= torque efficiency 
= effective bulk modulus 
= moment of inertia of payload 
= moment of inertia of 
= mode shape coefficient 
= mode shapew coefficient 
= computed mode shape coefficient 
= computed mode shape coefficient 
= to evaluate the coefficient and be used in fa 

= to evaluate the coefficient and be used in maa 

= to evaluate the coefficient and be used in mab or mba 

= to evaluate the coefficient and be used in maa 

= to evaluate the coefficient and be used in kb 

= coefficient of large motion acceleration in large motion dynamics 
= coefficient of small motion acceleration in large motion dynamics 
= coefficient of large motion acceleration in small motion dynamics 
= coefficient of smalle motion acceleration in small motion dynamics 
= coefficient of stiffness in small motion dynamics 
= right-hand-side of large motion dynamics 
= right-hand-side of small motion dynamics 
= large motion angle and thetaold represents previous one in 
computation, thetanew represents the new one. 
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* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

♦ 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 


thetad = velocity of large motion with the same expression as theta 
thetadd = acceleration of large motion with tha same expression as theta 

v = small motion and void represents previous one in the computation 

and the vnew represents the new one. 

vd = velocity of small motion and has the same expression as in v. 

vdd = acceleration of small motion and has the same expression as in v 

pi = load hydraulic pressure drop 

Ps = hydraulic supply pressure 

current = the current input the electrohydraulic actuactor 
psi = total angle is the sum of large motion (theta) and small motion (v) 

psiv = velocity of total angle 

psia = acceleration of total angle 

psir = desired total angle 

psivr = the first time derivative of psir 

psiar = the second time derivative of psir 

psijk = the third time derivative of psir 

Function expalnation 

constl = to perform the computation of constant in the program 
timinvarm = to perform the computation of time invariant terms 

ic = to perform the computation of initial condition 

tracking = to give the desired trajectory 

sequenm = to perform the sequential integration method 

tippva = to perform the transformation of total expression 

ismco = to compute the current using sliding control with first-order 

sliding condition 

ismcst = to compute the current using the straight sliding control 

ismcva = to compute the current using the versatile sliding control 

hydralm = to perform the dynamics of electrohydraulic actuactor 

plotter = to perform the plot 




% sliding control with first-order sliding condition algorithm 
clear,clg 

tra=input(’Enter desired trajectory 1 or 2 ’) 
lambda=input(’please enter lambda ’) 

h=0.001 ;FT=10;t=0:h:FT; 

ml=0;speed = 10;R=[300;30;l];G=[1000;0;0]; 

Bhat=[0;0;l.le4];DB=295;fhat=[0;0;-5.7732e4];Df=6.7117e4; 

Ahat=[0 1 0;0 0 1;0 0 0j;DA=0; 

[L,rho,A,E,I,Dv,g,K,Dm,Vt,Ps,Ctm,eff,betae,Ip,Ir,betal ,beta2,Cl ,C2,B1 ,B2] 
=constl(ml); 
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% computation for the time invariant terms 

[nlm,nlnlm,lxnlm,swwsm,klei,mab,mba,mbb]=timinvarm(rho,A,E,I,betal,beta2. 

,Cl,C2,Bl,B2,Dv,Ip,L,ml); 

% Initial Conditions 

[thetaold,thetadold,thetaddold,vold,vdold,vddold,torqueactold,pl,current(l).... 
,psi(l),psiv(l),psia(l)]=ic(Dm,eff,Ctm,Ps,K,L); 

% desired trajectory 
[psir,psivr,psiar,psijk]=tracking(t,tra); 

for n=1: (length(t)-1) 

n % sequential integration method 

[thetanew,thetadnew,thetaddnew,vnew,vdnew,vddnew]=sequenm(thetaold,thetadold... 

,thetaddold,vold,vdold,vddold,mab,mba,mbb,swwsm,nlnlm,lxnlm,nlm,klei,Ip. 

,Ir,Dv,ml,g,h,L,A,rho,torqueactold); 

[psi(n+ l),psiv(n + l),psia(n+1)]=tippva(thetanew,thetadnew,thetaddnew,vnew. 

, vdnew, vddnew, L); 

% sliding mode control with first-order sliding condition 

[current(n+l),s(n+ l),intgerror,phi(n+1)] =ismco(Ahat,Bhat,fhat,DA,DB,Df,R,G.... 

,h,lambda,psi(n+l),psiv(n+l),psia(n+l),psir(n+l),psivr(n+l),psiar(n + l). 

,psijk(n+1),speed,intgerror); 

% electriohydraulic actu actor 

[plnew,torqueactnew] =hydralm(h,current(n + l),thetadnew,pl,K,Dm,Vt,Ps,Ctm. 

,eff,betae); 
pi = plnew; 

% 

thetaold=thetanew;thetadold=thetadnew;thetaddold=thetaddnew; 
void=vnew; vdold=vdnew; vddold = vddnew; 
torqueactold = torqueactnew; 
end 

% plotter - > 

[ps,pphi,pciirrent,ppsi]=plotter(s,phi,current,psi) 

******************•***********•****************************•****,,***********•* 


% Straight Sliding Mode Control Algorithm 
clear, clg 

tra=input(’Enter desired trajectory 1 or 2 ’) 


89 







w=input(’please enter w ’) 
zeta=input(’please enter zeta ’) 


h=0.0002;FT=1.2;t=0:h:FT;ml=0; 
R=[300;30; l];G=[1000;0;0];speed=10; 


% Uncertainties 

Bhat=[0;0; 1. le4];DB=295; 

fhat=[0;0;-5.7732e4];Df=6.7114e4; 

Ahat=[0 1 0;0 0 1;0 0 Q];DA=0; 

% Define constants 

[L,rho,A,E,I,Dv,g,K,Dm,Vt,Ps,Ctm,eff,betae,Ip,Ir,betal,beta2,Cl,C2,Bl,B2]. 

=const 1 (ml); 

% computation for some time invariant terms 

[nlm,nlnlm,lxnlm,swwsm,klei,mab,mba,mbb]=timinvarm(rho,A,E,I,betal,beta2 

,Cl,C2,Bl,B2,Dv,Ip,L,ml); 

% Initial conditions 

[thetaold, thetadold, thetaddold, void, vdold, vddold, torqueactold ,pl. 

,current(l) ,psi(l) ,psiv( 1) ,psia(l)] =ic(Dm,eff,Ctm,Ps,K,L); 

% desired trajectory 
[psir,psivr,psiar,psijk] =tracking(t,tra); 

% Define initial values used in straight and intgral error sliding control 

[intgerr,intgs,s(l),phi(l)]=smciv(R,G,psi(l),psiv(l),psia(l),psir,psivr. 

,psiar,h,w,zeta); 

for n=1: (length(t)-1) 
n 

[thetanew,thetadnew,thetaddnew,vnew,vdnew,vddnew]=sequenm(thetaold. 

.thetadold,thetaddold,void,vdold,vddold,mab,mba,mbb,swwsm,nlnlm,lxnlm,nlm. 
,klei,Ip,Ir,Dv,ml,g,h,L,A,rho,torqueactold); 

[psi(n+ l),psiv(n+l),psia(n+1)] =tippva(thetanew,thetadnew,thetaddnew. 

,vnew,vdnew,vddnew,L); 

% straight sliding mode control 

[current(n+ l),s(n+ l),intgs,intgerr,phi(n+1)] =ismcst(Ahat,Bhat,fhat,DA,DB. 

,Df,R,G, h,w,psi(n +1) ,psiv(n +1) ,psia(n +1) ,psir(n +1) ,psivr(n +1) ,psiar(n +1). 

,psijk(n +1), speed, intgs, in tgerr, zeta); 
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% electrohydraulic actuactor 

[plnew,torqueactnew]=hydralm(h,current(n+l),thetadnew,pl,K,Dm,Vt,Ps,Ctm. 

,eff,betae); 

pl=plnew; 

% 

thetaold=thetanew;thetadold=thetadnew;thetaddold=thetaddnew; 
void=vnew; vdold —vdnew;vddold=vddnew; 
torqueactold=torqueactnew; 
end 

% plotter -- > 

[ps,pphi,pcurrent,ppsi]=plotter(s,phi,current,psi) 

****************************************************************************** 


% Versatile Sliding Mode Control Algorithm 
clear,clg 

tra=input(’Enter desired trajectory 1 or 2 ’) 

w=input(’please enter w ’) 
zo=input(’please enter z ’) 
zeta=input(’please enter zeta ’) 

h-0.0002;FT= 1.2t=0:h:FT;ml=0; 
speed = 10;R=[300;30;l];G=[1000;0;0]; 

% Uncertainties 

Bhat=[0;0;l.le4];DB=295; 

fhat=[0;0;-5.7732e4];Df=6.7117e4; 

Ahat=[0 1 0;0 0 1;0 0 0];DA=0; 
s=0;intgerr=0; 

% Define constants 

[L,rho, A,E,I,Dv,g,K,Dm,Vt,Ps,Ctm,eff,betae,Ip,Ir,betal ,beta2,Cl ,C2,B1 ,B2]. 

=const 1 (ml); 

% computation for some time invariant terms 

[nlm,nlnlm,lxnlm,swwsm,klei,mab,mba,mbb]=timinvarm(rho,A,E,I,betal,beta2 
,C1 ,C2,B1 ,B2,Dv,Ip,L,ml); 

% Initial conditions 

[thetaold,thetadold,thetaddold,void,vdold.vddold,torqueactold,pi,current(l). 

,psi(l),psiv(l),psia(l)]=ic(Dm,eff,Ctm,Ps,K,L); 


% Desired trajectory 
[psir,psivr,psiar,psijkj=tracking(t,tra); 
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for n=l:(length(t)-l) 
n % sequential integration method 

[thetanew, thetadnew, thetaddnew, vnew, vdnew, vddnew]=sequenm(thetaold. 

,thetadold,thetaddold,vold,vdold,vddold,mab,mba,mbb,swwsm,nlnlm,lxnlm. 

,nlm,klei,Ip,Ir,Dv,ml,g,h,L,A,rho,torqueactold); 

[psi(n+l),psiv(n+l),psia(n+1)] =tippva(thetanew, thetadnew, thetaddnew. 

,vnew,vdnew,vddnew,L); 

% Versatile sliding mode control 

[current(n+l),sdot(n),s,intgerr,phi(n)]=ismcva(Ahat,Bhat,fhat,DA,DB,Df,R,G.... 

,psi(n+ l),psiv(n+ l),psia(n + l),psir(n+l),psivr(n+ l),psiar(n+l),psijk(n+1). 

,speed,s,intgerr,h,w,zo,zeta); 

% electrohydraulic actuator 

[plnew,torqueactnew] =hydralm(h,current(n+1),thetadnew,pi,K,Dm, Vt,Ps. 

,Ctm,eff,betae); 

pl=plnew; 

% 

thetaold=thetanew; thetadold=thetadnew; thetaddold=thetaddnew; 
void=vnewjvdold=vdnew;vddold=vddnew; 
torqueactold=torqueactnew; 
end 

% plotter - > 

[psdot,pphi,pcurrent,ppsi] =plotter(sdot, phi, current, psi) 

m***************************************************************************** 


function [L,rho,A,E,I,Dv,g,K,Dm,Vt,Ps,Ctm,eff,betae,Ip,Ir,betal,beta2,Cl 
,C2,B1,B2] = constl(ml) 


rho = 7861.05; A = 6.17795e-04; 

L = 0.9985;g = [0;0;-9.8066]; 

E = 2.0ell;I = 4.065e-10; 

Dv= [0 0 0;0 0 0; 1 0 0]; % first derivative of D w.r.t. v 

Iyy = (1.3653333e-3)*ml; 

Ixx= ((ml/0.4233) A 2)*(7.5e-06)*ml; 

Ip = [ml 0 0;0 Ixx 0;0 0 Iyy]; % Ip = moment of inertia of payload 
mr = 9.00011451;Ir = [mr 0 0;0 0.02746713 0;0 0 0.02746713]; 

K = 2.402963e-09; % unit = ( m A 4/(sec-mA-sqrt(N» ) 

Dm = 6.227le-05; % unit = ( m*3/rad ) 

Vt = 3.05127e-04; % unit = ( m A 3 ) 

Ps *= 1.37888e07; % unit = ( N/m A 2 ) 

Ctm = 3.7064772e-13; % unit = ( m A 5/N-sec ) 
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betae = 6.90e08; % unit = ( N/m A 2 ) 

eff = 0.9; 

betalL = 1.875104069;beta2L = 4.694091133; 
betal = betalL/L;beta2 = beta2L/L; 

Cll = sin(betalL) + sinh(betalL); 

C12 = cos(betalL) + cosh(betalL); 

Cl = C11/C12; 

C21 = sin(beta2L) + sinh(beta2L); 

C22 = cos(beta2L) + cosh(beta2L); 

C2 « C21/C22; 

D = 4*Cl*beta2 - 4*C2*betal; 

B1 = 2*beta2/D;B2 = -2*betal/D; 

***********************^*************************%**************************** 


function [alfa,delta,gamaO,gamal,gama2,gama3,gama4,gama5] = const2(h) 

alfa = 0.25;delta = 0.5; 

gamaO = l/(alfa*h A 2); 
gamal = delta/(alfa*h); 
gama2 = l/(alfa*h); 
gama3 = (l/(2*alfa))-l; 
gama4 = (delta/alfa)-l; 
gama5 = (h/2)*((delta/aIfa)-2); 

****************************************************************************** 


function [alfa,delta,gamaO,gamal,gama2,gama3,gama4,gama5,gama6,gama7] =const3(h) 

alfa = 0.25;delta = 0.5; 

gamaO = l/(aIfa*lT2); 
gamal = delta/(alfa*h); 
gama2 - l/(alfa*h); 
gama3 = (l/(2*alfa))-l; 
gama4 = (delta/alfa)-l; 
gama5 = (h/2)*((delta/alfa)-2); 
gama6 = h*(l-delta); 
gama7 = delta *h; 

****************************************************************************** 


function [nlm,nlnlm,lxnlm,swwsm,klei,mmab,mmba,mmbb]=timinvarm(rho,A,E,I 
,betal,beta2,Cl,C2,Bl,B2,Dv,Ip,L,ml) 

swwsm = simps(’swws\0,0,0,0,0,0)*rho*A; 
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nlm = $imps(’nr,betal,beta2,Cl,C2,Bl*B2)*rho*A; 

klei = simps(’kl 1 ’ ,betal ,beta2,C 1 ,C2,B1 ,B2)*E*I; 
nlnlm = simps(’nlnr,betal,beta2,Cl,C2,Bl»B2)*rho*A; 
lxnlm = simpsClxnl’, betal,beta2,Cl,C2,Bl>B2)*rho*A; 

[D,Dd,w,wd,wdd,Ad,Add] = usualm(-pi/3,0,7,L); 
mmba = mba(D,w,wd,Dv, Ip, lxnlm); 
mmbb = mbb(w,Dv,Ip,nlnlm,ml); 
mmab = mab(D,w,wd,Dv,Ip,lxnlm,ml,L); 

m++mm*+*+##m+++**+m*+mm+m***++***m+**+*+++*+m*++***********+********#***m***** 


function [thetal,thetadl,thetaddl,vl,vdl,vddl,torqueactl,pl,currentl,psil 
,psivl,psial]=ic(Dm,eff,Ctm,Ps,K,L) 

thetal =0;thetadl =0;thetaddl =0; 

vl =-0.031 ll;%vl =-0.061 ll;%vl =-0.146064149; %vl =-0.23155; 
vdl=0;vddl=0; 


torqueactl = 44.5549095;pl=torqueactl/(Dm*eff); 

Q=Ctm*pl;I0=Q/(K*sqrt(Ps-pl)) ;Imax =10; 
ifirac=0.5*(Imax-I0);currentl =I0+ifrac; 

m**m*********** *************************************************************** 


% Trajectory design function 
% 

function [psid,psivd,psiad,psijerk]=tracking(t,n) 


if (n = = l) 

psid =ones(t);psivd=zeros(t);psiad=zeros(t);psijerk=zeros(t); 
elseif(n==2) 

psijerk=0;psiad=0; 
for n = l:length(t) 
if (t(n)< 1) 

psid(n)=t(n);psivd(n) = 1; 
elseif (1 < =t(n) & t(n)<2) 
psid(n) = 1 ;psivd(n) =0; 
end 

if (2< =t(n) & t(n) < 2.4) 
psid(n)=3-t(n) ;psivd(n)=-1; 
elseif (t(n)> =2.4) 

psid(n)=0.6;psivd(n)=0; 

end 

end 

end 

+**+***mm*****+m****+*+*m**+++*++**+*++**+*+»**mm+*+*m+***+m***++************* 
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% 

% Sequential integration method is used to do numerical iteration. 

% 

function [thetanew,thetadnew,thetaddnew,vnew,vdnew,vddnew] = sequen(theta. 

,thetad,thetadd,v,vd,vdd,mab,mba,mbb,swwsm,nlnlm,lxnlm,nlm,klei,Ip,Ir,Dv. 

,ml,g,h,L,A,rho,torque) 

[alfa,delta,gama0,gamal,gama2,gama3,gama4,gama5] = const2(h); 

% following steps are solving for small-motion(vnew, vdnew, & vddnew) 

% (1) original coefficients computation (large and small-motion ) at time t 
[D,Dd,w,wd,wdd,Ad,Add] = usualm(theta,v,vd,L); 

Kb *= kb(thetad,klei,nlnlm); 

Fb = fb(thetad,w,wd,D,Dd,nlm,Dv,Ip,g,ml); 

Maa = maa(v,wd,Ad,D,Ip,Ir,swwsm,nlnlm,A,rho); 

Fa = fa(v,vd,thetad, wd, wdd, Ad, Add,D,Dd,Ip,Ir,nlm,nlnlm,g,ml,torque); 

% (2) modify small-motion coefficients ( at time t) to avoid secular terms 
Mbs = mbb - mba*(l/Maa)*mab; 

Fbs = Fb - mba*(l/Maa)*Fa; 
omega = (l/(2*pi))*sqrt (Kb/Mbs); 

C = 0.2*(sqrt (Kb/Mbs)); 

% (3) fast (small) motion solver ~> vnew, vdnew, and vddnew 
Kbi = Kb + gamaO*(Mbs) 4- gamal*C; 

Fbi=Fbs+Mbs^gamaO* v+gama2 *' vd+gama3 * vdd)+C * (gama 1 *v+ gama4 *’ vd+gama5 *■vdd); 
vnew = (inv(Kbi))*Fbi; 

vdnew=vd+(l-delta)*h*vdd+gamal , '(vnew-v)-((delta/alfa)*vd)-gama3 ,, ‘delta*h*vdd; 
vddnew = gamaO*(vnew-v) - gama2*vd - gama3*vdd; 

% using predictor and corrector scheme to find large-motion 
% (1) using predictor and new etas’ to evaluate large-motion coefficients 
thetap = theta + thetad*h + (0.5-alf?.)*(h*h)*(thetadd); 
thetadp = thetad + (1-delta) *h*theiadd; 

[D,Dd,w,wd,wdd,Ad,AddJ = usualm(thetap,vnew,vdnew,L); 

Maa = maa(vnew,wd,Ad,D,I^,Ir,swwsm,nlnlm,A,rho); 

Fa = fa(vnew,vdnew,thetadp, wd, wdd, Ad, Add, D,Dd,Ip,Ir,nlm,nlnlm,g,ml, torque); 

% (2) solve for thetaddnew 
Fae = Fa - mab“vddnew; 
thetaddnew = Fae/(Maa); 

% (3) find correctors using thetaddnew 
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thetac = alfa*(h*h)*thetaddnew; 
thetadc = delta*h*thetaddnew; 

% (4) solve for thetanew and thetadnew 
thetadnew = thetadp + thetadc; 
thetanew = thetap + thetac; 

****************************************************************************** 


function [x,xr,xrd]=tippva(theta,thetad,thetadd,v,vd,vdd,psid,psivd,psiad. 

,psij,L) 

psi = theta + (v/L); 
psiv = thetad + (vd/L); 
psia = thetadd + (vdd/L); 
x=[psi;psiv;psia]; 
xr—[psid;psivd;psiad]; 
xrd = [psivd;psiad;psij]; 

****************************************************************************** 


function [currents,intgerror,phi]=ismco(Ahat,Bhat,fhat,DA,DB,Df,R,G,h,lambda. 
,psi,psiv,psia,psir,psivr,psiar,psijk,speed,intgerror); 

x=[psi;psiv;psia] ;xr=[psir;psivr;psiar] ;xrd=[psivr;psiar;psijk]; 

error = x - xr; 

intgerror=intgerror+error*h; 

RtB=inv(R’*Bhat);RtBn=norm(RtB,2); 
s=R’ *error+G’ *intgerror; 

RtBsgn=RtB*sign(s);RtBsn=norm(RtBsgn ,2); 

% norminal value 

Ihat = R’ *(xrd-fhat-Ahat*x)-G’ *error; 

Ihatn= norm(Ihat,2); 

boxl =RtB*Ihat; 
box 1 n=norm (box 1,2); 

% k 

Rn=norm(R,2); 

k=((Rn)*(Df + DB*box 1 n)+speed)/( 1 -(DB*'RtBsn)); 
phi=(k/lambda) *( 1 + DB^RtBn); 


std — s/phi; 
if( abs(s) > =phi ) 
if( std > — 0 ) 
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sat=l; 

else 

sat=-1; 
end 

elseif( abs(s) < phi) 
sat=std; 
end 

current=(RtB) *(Ihat-k*sat); 
if( current > = 10) 
current =10; 
end 

if( current < = 10) 
current=-10; 
end 


function [current,s,intgs,intgerr,phi]=ismcst(Ahat,Bhat,fhat,DA,DB,Df,R,G 
,h,w,psi,psiv,psia,psir,psivr,psiar,psijk,speed,intgs,intgerr,zeta); 

x = [psi;psiv;psia] ;xr=[psir;psivr;psiar] ;xrd=[psivr;psiar;psijk]; 

error=x-xr; s=R ’ ♦error; 

intgs=intgs+s*h; 

RtB=inv(R’ ♦Bhat) ;RtBn=norm(RtB,2); 

RtBsgn=RtB*sign(s);RtBsn=norm(RtBsgn,2); 

% nominal value 

Ihat=R’*(xrd-fhat-Ahat*x)-(w^2)*intgs; 

Ihatn=nom(Ihat,2); 

boxl=RtB*Ihat; 
box 1 n=norm (box 1,2); 

% k 

Rn=norm(R,2); 

k=((Rn)*(Df+DB*box 1 n)+speed)/( 1 -(DB*RtBsn)); 

% tuning parameters 

lambda=(sqrt(l +DB*RtBn))*w; 

phi=(k/(2 ,> zeta*lambda)) , *(DB*RtBn); 


std — s/phi; 
if( abs(s)> =phi ) 
if( std> =0 ) 
sat=l; 
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else 

sat=-1; 
end 

elseif( abs(s)<phi) 
sat=std; 
end 

% (3) current 
current=RtB*(Ihat-k*sat); 
if( current > = 10) 
current=10; 
end 

if( current < =-10) 
current=-10; 
end 

*************************************************t**************************** 


function [current,sdot,s,intgerr,phi] =ismcva(Ahat,Bhat,fhat,DA,DB,Df,R,G 
,psi,psiv,psia,psir,psivr,psiar,psijk,speed,s,intgerr,h,w,zo,zeta); 

x=[psi;psiv;psia] ;xr=[psir;psivr;psiar] ;xrd=[psivr ;psiar;psijk]; 
error=x-xr; 

intgerr=intgerr+error*h; 

sdot=R’ *error+G’ *intgerr-zo*s; 

s=s+sdot*h; 

RtB=in v(R ’ *Bhat); RtBn=norm(RtB, 2); 

RtBsgn =RtB*sign(sdot);RtBsn=norm(RtBsgn,2); 

% norminal value 

Ihat=R’ *(xrd-fhat-Ahat*x)-(w A 2)*s-G’ *error+zo*sdot; 

Ihatn=norm(Ihat,2); 

boxl =RtB*Ihat; 
boxln=norm(boxl ,2); 

% k 

Rn=norm(R,2); 

k=((Rn)*(Df+DB*boxln)+speed)/(l-(DB*RtBsn)); 

% tuning parameters 

lambda=(sqrt(l +DB*RtBn))*w; 

phi=(k*(DB*RtBn))/((2 *zeta*lambda)+(DB*RtB*zo)); 

% sat(sdot) FUNCTION 
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std=sdot/phi; 
if( abs(sdot) > =phi) 
if( std> =0) 
sat=l; 
else 

sat=-l; 

end 

elseif( abs(sdot)<phi) 
sat=std; 
end 

% CURRENT 
current=RtB*(Ihat-k*sat); 
if( current > = 10) 
current =10; 
end 

if( current < =-10) 
current=-10; 
end 

** 4 ^****** 4 ^***** 4 ^*********************** 4 ^********************************** 


function [Plnew,torque] =hydralm(h,current,thetad,PI,K,Dm,Vt,Ps,Ctm,eff,betae) 
% 

% A fourth-order Runge-Kutta method is used to obtain the numerical solution 
% In this case, the only variable will be PI 
% 

Pv = ’Ps - (Pl+x)’; 

Pld = ’(4*betae/Vt) * (K*current*(sqrt(eval(Pv))) - Dm*thetad - Ctm*(Pl-t-x))’; 

% To evaluate kl, let x = 0 
x = 0; 

if (eval(Pv)) < 0 
Pv = ’O’; 
end 

kl = h*eval(Pld); 

% To evaluate k2, let x = kl/2 
x = kl/2; 
if (eval(Pv)) < 0 
Pv = ’O’; 
end 

k2 = h*eval(Pld); 

% To evaluate k3, let x = k2/2 
x = k2/2; 
if (eval(Pv)) < 0 
Pv = ’O’; 
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end 

k3 = h*eval(Pld); 

% To evaluate k4, let x = k3 
x = k3; 

if (eval(Pv)) < 0 
Pv = ’O’; 
end 

k4 = h*eval(Pld); 

Plnew = PI + (l/6)*(kl 4- 2*k2 + 2*k3 + k4); 
ifPlnew > Ps 
Plnew = Ps; 
end 

ifPlnew < -Ps 
Plnew = -Ps; 
end 

torque = eff*Plnew*Dm; 

****************************************************************************** 


function [ps,pphi,pcurrent,ppsi] =plotter(s,phi,current,psi,t) 
subplot(221) ,plot(t,s),xlabel(’Time(sec)’) ,ylabel(’S ’) 
subplot(221),plot(t,phi),xlabel(’Time(sec)’),ylabel(’Phi’) 
subplot(221),plot(t,current),xlabel(’Time(sec)’),ylabel(’Current (mA)’) 

subplot(221) ,plot(t,psi) ,xlabel(’Time(sec) ’) ,ylabel(’Psi (Rad) ’) 

****************************************************************************** 


function [D,Dd,w,wd,wdd,Ad,Add] =usualm(tha,v,vd,L) 

D=[l 0 0;0 1 0;v 0 l];Dd=[0 0 0;0 0 0;vd 0 0]; 
w=[l 0 0;L*cos(tha) cos(tha) -sin(tha);L*sin(tha) sin(tha) cos(tha)]; 
wd=(0 0 0;-L*sin(tha) -sin(tha) -cos(tha);L*cos(tha) cos(tha) -sin(tha)]; 
wdd=[0 0 0;-L*cos(tha) -cos(tha) sin(tha);-L*sin(tha) -sin(tha) -cos(tha)]; 

Ad=[0 0 0;0 -sin(tha) -cos(tha);0 cos(tha) -sin(tha)); 

Add=[0 0 0;0 -cos(tha) sin(tha);0 -sin(tha) -cos(tha)]; 
****************************************************************************** 


function ans = fa(v,vd,thetad,wd,wdd,Ad,Add,D,Dd,Ip,Ir,nlm,nlnlm,g,ml,torque) 

% 

% evaluation of right hand side of the large-motion equation 
% 

% (1) -2*thetad*etad’*{integration [N(x)’*W,’*W,*N(x)] dm}*eta 
% where {integration [N(x)’*W,’*W,*N(x)]dm} = [nwwnl nwwnl2;nwwnl2 nwwn2] 

x = -2*thetad*vd*nlnlm*v; 
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% (2.1) (integration[s’] dm}*W,’*g 

% where {integration[s’] dm} = s <-- 1st moment of inertia of link 
s = [4.856519 -2.428258693 0]; 
yl = s*wd’*g; 

% (2.2) eta’*{integration[N’] dm}*W,’*g 
% where {integration[N’] dm} = n 

n = [0 0 nlm]; 

y2 = v*n*wd’*g; 

% (2.3) {integrationtSp’ldmJ^D’^W/^g 

% where {integration[Sp’] dm} = Sp <- 1st moment of inertia of payload 
Sp = [ml 0 0]; 
y3 = Sp*D , *wd”*‘g; 
y = yl+y2+y3; 

% (3) evaluation "trace" terms of right hand side of large-motion equation 

% trace [2 * wd *D *Ip*Dd ’ * wd ’ ]*thetad 

zll = -2*thetad , ' , wd*D ,, ■Ip*Dd , *wd , ; 

z = trace(zll); 

ans = x 4-y 4-z 4-torque; 

****************************************************************************** 


function ans = fb(thetad,w,wd,D,Dd,nlm,Dv,Ip,g,ml) 

% 

% evaluante the right hand side of the small-motion equation 
% 

% (1) {Integration [N’] dm}*W’*g 
n = [0 0 nlm]; 
u = n*w'*g; 
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% ( 2 ) 

% [ {integration [Sp’] dm}*Dv’*w’*g {integration [Sp’] dm}*Dphi’*w’*g ] 

% where {integration [Sp’] dm} = Sp 

Sp = [ml 0 0]; 

x = Sp*Dv’*w’*g; 

% (3) 

% evaluate all the "trace" terms of right hand side of small-motion equation 
% 

% (3.1) trace [wdd*D*Ip*Dv’*w’*(thetad A 2))+2*wd*Dd*Ip*Dv’ *w ’] 

yll = wd*Dv*Ip*D’ *wd’ *(thetad A 2); 
y22 = 2*wd*Dv*Ip*Dd’*w’*(thetad); 
y = trace(yll+y22); 

ans = u+x+y; % sum of the right hand side of small-motion equation 

****************************************************************************** 


function ans = maafv.wd.Ad.D.IpJr.swwsm.nlnlm.A.rho) 

% % 

% coefficient computation of "thetadd" of large-motion equation % 

% % 

% (1) {integration [s’*wd , *wd*s] dm} = swwsm 
x = swwsm; 

% (2) eta’*{integration [N’*wd’*wd*N] dm}*eta 
y = v’*nlnlm*v; % 2nd-order value should be small 

% The following is used to evaluate all the "trace" terms in the coefll 
% (2.1) Trace [ W^D^Ip^D^W,’ ] 
ul = trace(wd *D *Ip*D ’ *wd ’) ; 

% (2.2) Trace [ A,*Ir*A, ] 
u2 = trace(Ad*Ir*Ad’); 
z = ul+u2; 

ans = x-fy+z; 


function ans = mab(D,w,wd,Dv,Ip,lxnlm,ml,L) 

% 

% coefficient computation of "etadd" of large-motion equation 
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% 

% (1) {integration [s’*wd’*w’ ,, N] dm} = [swwnl swwn2] 
x = lxnlm; 

% (2) [trace (W,*D*Ip*(Dv)’*W’) trace (W,*D*Ip*(Dphi)’*W’)] 
y = trace(wd*D*Ip*Dv’*w’); 

ans = x+y; 

******************************************** i********************************** 

function ans = mba(D,w,wd,Dv,Ip,lxnlm) 

% 

% coefficient computation of "etadd" of large-motion equation 
% 

% (1) {integration [s’*wd”*‘w*N] dm} = [swwnl swwn2] 
x = lxnlm; 

% (2) [trace (W,*D*Ip*(Dv)’*W’) trace (W,*D*Ip*(Dphi)’*W’)] 
y = trace(wd*D*Ip*Dv’*w’); 

ans = x+y; 

♦a**************************************************************************** 

function ans = mbb(w,Dv,Ip,nlnlm,ml) 

% 

% coefficient computation of "etadd" of small-motion equation 
% 

% (1) {integration [N , *w*w*N] dm} = [nwwnl nwwnl2;nwwnl2 nwwn2] 
x = nlnlm; 

% (2) [ Tr(W*Dv*Ip*Dv’*W’) 0;0 Tr(W*(Dphi)*Ip*(Dphi)’*W’) ] 
y = trace(w*Dv*Ip*Dv’*w’); 

ans = x+y; 

********** 4 ^****************************************************************** 

function ans = kb(thetad,klei,nlnlm) 

% 

% coefficient computation of "eta" of small-motion equation 
% 

% (1) ((E*I)/2)*{integartion [sigma+sigma’] dx} = x 
x = klei; 

% (2) {integration [N’^wd’^wd^N] dm}*eta*'(thetad A 2) 
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% where {integration [N’*wd’*‘wd*N] dm} = [nwwnl nwwnl2;nwwnl2 nwwn2] 
y = -(thetad A 2)*nlnlm; 


ans=x+y; 

mm**************************************************************************** 


function ans = nl(x,betal,beta2,Cl,C2,Bl»B2) 

% 

% natural-mode shape functions Nl(x) evaluation 
% 

xil = Cl*(cos(betal*x)+cosh(betal*x)) + (sin(betal*x)+sinh(betal*x)); 
xi2 = C2*(cos(beta2*x)+cosh(beta2*x)) + (sin(beta2*x) + sinh(beta2*x)); 
ans = Bl*xil + B2*xi2; 

nmt+m************************************************************************* 


function ans = nlnl(x,betal,beta2,Cl,C2,Bl,B2) 

% 

xil = Cl*(cos(betal*x)+cosh(betal*x)) + (sin (beta l*x)+sinh (beta l*x)); 
xi2 = C2*(cos(beta2*x)+cosh(beta2*x)) + (sin(beta2*x)+sinh(beta2*x)); 

N1 = Bl*xil + B2*xi2; 

% ’nwwnl* is to evaluate N1 A 2 and to be used in "coefll.m M 
% An evaluation of [ N1 A 2 N1*N2 ; N1*N2 N2 A 2 ] which comes from multiplic- 
% ation of [ N(x)’*W,’*W,*N(x) ] 
ans = N1. A 2; 

****************************************************************************** 


function ans = lxni(x,beial,beta2,Cl,C2,Bl,B2) 

xil = Cl*(cos(betal*x)+cosh(betal*x)) + (sin(betal*x)+sinh(betal*x)); 

xi2 = C2*(cos(beta2*x)+cosh(beta2*x)) + (sin(beta2*x)+sinh(beta2*x)); 

N1 = Bl*xil + B2*xi2; 
ans = 0.9985 *N1 + x.*Nl; 

****************************************************************************** 


function y = swws(x,betal,beta2,Cl,C2,Bl,B2) 

% 

% To be used in "coefll.m" and evaluate integral [S’*W,’*W,*S]dx 
% After multiplication [ S’*W,’*W,*S ] = (L+x) A 2 
% where W represents 2-D homogeneous transform matrix 
% W, represents first derivative w.r.t. theta 
% S represents the local position along the link [1 x 0]’ 
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% L length of the link = 0.9985 m 
% 

y = (0.9985+x). A 2; 

mm*+***+m+++*m#m++m+*+m+m*+*#+++**********+*+********************************* 


function ans = kll(x,betal,beta2,Cl,C2,Bl>B2) 

xi 1 dd=(beta 1 *beta 1) * (C1 *(-cos(beta 1 *x) + cosh (beta 1 * x))+(-sin(beta 1 * x) 
+sinh(betal*x))); 

xi2dd = (beta2 *beta2) *(C2 *(-cos(beta2 *x)+cosh(beta2 *x))+(-sin(beta2 *x) 
+sinh(beta2*x))); 

% Note: Nldd = second derivative of Nl(x) w.r.t. x 
% N2dd = second derivative of N2(x) w.r.t. x 
Nldd = Bl*xildd + B2*xi2dd; 
ans = Nldd.^2; 
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APPENDIX B 


FORTRAN IMPLEMENTATION PROGRAMS 


C First-Order Sliding Condition + An Integral Error Control 
SINCLUDE: ’ATLDEFS.FOR’ 

SINCLUDE: ’ ATLERRS .FOR ’ 

PARAMETER (NUM = 1750) 

INTEGERS ADGAINS( 16), ADCHAN( 16), ICONFIG( 16) ,B ASE ADR, CDE VID, 
+CDEVF,LG,SCAN,DA VAL1,DAVAL2,Y,STATUS 
INTEGERS FH,FM,FS,FSS 

REAL BHAT,FHAT,L,PI,K,R(3),G,IHAT,INTGERR,LAMBDA 

REAL PSID(NUM),PSIVD(NUM),PSUERK(NUM),S(NUM),CURRENT(NUM) 

REAL PSI(0:NUM),PHI(NUM) 

COMMON/CONFIG/ICONFIG 

EQUIVALENCE 

(ICONFIG(KCBASEADR), B ASEADR), (ICONFIG(KCDEVID) ,CDEVID) 

+, (ICONFIG(KCDEVFLAGS),CDEVFLG), (IC0NFIG(KCSCAN) ,SC AN),(ICONFIG(KCC 
+H ANNELS), CHAN) 

C 

DATA R,G/300.0,30.0,1.0,3000.0/,SPEED,L/10.0,0.9985/ 

DATA SMALFIL,INTGERR/2*0.0/,RATE,FREQ/2000.0,200.0/ 

DATA PSUERK/NUM*0.0/,PSrVOD,D20LD,VOLD/3*0.0/ 

C UNCERTAINIES 

DATA BHAT,DB/1.1E4,295.0/,FHAT,DF/-5.7732E4,6.7114E4/ 

C 

PRINT VPlease Enter Lambda for S dynamics, and Payload’ 

READ *,LAMBDA,PAYLOAD 

PRINT *,’Please Enter Desired Trajectory "l",or "2" ’ 

READ *,IANS 
C DESIRED TRAJECTORY 
IF(IANS .EQ. 1) THEN 
DO 100 M = 1,NUM 
TIME=REAL(M)/FREQ 
IF (TIME .LE. 1.0) THEN 
PSID(M) =TIME 
PSIVD(M) = 1.0 

ELSE IF((TIME .Gc. 1.0) .AND. (TIME .LE. 2.0)) THEN 
PSID(M) = 1.0 
PSIVD(M)=0.0 
END IF 

IF((TIME .GT. 2.0) .AND. (TIME .LE. 2.4)) THEN 
PSID(M)=3.0-TIME 
PSIVD(M) = -1.0 
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ELSE IF(TIME .GE. 2.4) THEN 
PSID(M)=0.6 
PSIVD(M)=0.0 
END IF 

100 CONTINUE 

ELSE IF(IANS .EQ. 2) THEN 
DO 500 M = 1,NUM 
PSID(M) = 1.0 
PSIVD(M)=0.0 
PSUERK(M) =0.0 
500 CONTINUE 
END IF 
C 

M=0 

PI=4*ATAN(1.0) 

PAR=PI/1054.0 
R1B = 1.0/(R(3)*BHAT) 

RlNOM=SQRT(R(l)**2+R(2)**2+R(3)**2) 

OPEN(UNIT= 15,FILE=’C:\ATLAB\TMP\CTRL1 .M’,STATUS =’NEW’) 

C INITIALIZE DATA ACQUSITION BOARD 
STATUS = ALINIT0 
STATUS = ALSB(1) 

STATUS = ALSF(RATE) 

STATUS=ALRSET0 
STATUS=ALGC(ICONFIG) 

STATUS=ALDV(0,Y) 

C 

WRJTE(*,*)’ Are you ready (Ctrl-C) ? ’ 

WRITEO*,*)’ When you are ready to go, press "1" and "return" ’ 

READ MANSI 
IF(IANS 1 .EQ. 1) GOTO 7 
C BEGIN SIGNAL PROCESSING 
7 CALL GETTIM(IH,IM,IS,ISS) 

5 STATUS=ALAV(l,l,DAVAL1) 

STATUS = ALAV(2,1 ,DAVAL2) 

D1 =REAL(DAVALl-2048) 

D2=REAL(DAVAL2-2048) 

C FILTERED VERSION OF SMALL MOTION SIGNAL COMING OUT FROM STRAIN 
INDICATOR 

SMALFIL=0.94175*SMALFIL+0.029129*(D20LD+D2) 

C CALIBRATED EQUATION BETWEEN DISPLACEMENT(cm) AND SMALFIL 
SMALL=-0.3*(SMALFIL+15) 

SMANG=SMALL/(100.0*L) 

C LARGE MOTION CALIBRATED EQUATION 
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THETA=PAR*(D 1-30.0) 

C TOTAL ANGLE 

PSI(M) =THETA+SMANG 
C CONTROLLER DESIGN 
EF(M .EQ. 0) GOTO 1 
POSERR=PSI(M)-PSID(M) 

INTGERR=POSERR7FREQ + INTGERR 
C FILTERED VERSION OF TOTAL VELOCITY, PSIV(M) 

VNEW=(PSI(M)-PSI(M-1)) *FREQ 

PSIV=0.980198*PSIVOD+0.00990099*(VOLD+VNEW) 

VELERR=PSIV-PSIVD(M) 

S(M) =R( 1) "‘POSERR+R(2) *VELERR+G “"INTGERR 
C NORMINAL VALUE 

IHAT=R( 1) *(-VELERR)+R(3) "“(PS UERK(M)-FHAT)-G *POSERR 
C TO FIND K 

K=(RlNOM*(DF+DB*(RlB*IHAT))+SPEED)/(1.0-(DB*RlB)) 

PHI(M)=(K/LAMBDA)*(1.0+DB*RlB) 

STD=S(M)/PHI(M) 

EF(ABS(S(M)) .LT. PHI(M)) SAT=STD 
IF(ABS(S(M)) .GE. PHI(M)) THEN 
IF(STD .GE. 0.0) SAT=1.0 
DF(STD .LT. 0.0) SAT=-1.0 
ENDIF 

CURRENT(M)=R1B*((IHAT - K*SAT)/1000.0) 

IF (CURRENT(M) .GE. 10.0E-03) CURRENT(M) = 10.0E-03 
IF (CURRENT(M) .LE. -10.0E-03) CURRENT(M)=-10.0E-03 
VOLT=CURRENT(M) *500.0 

C DIGITAL TO ANALOG (D/A) CALIBRATED EQUATION 
Y=NINT (((V OLT+9.9868)/4.8766)*l000.0) 

STATUS =ALDV(0,Y) 

D20LD=D2 

VOLD=VNEW 

PSIVOD=PSIV 

C 

1 M=M+1 

IF ( M .LE. NUM ) GOTO 5 
WRn^WGame Over!’ 

CALL GETTIM(FH,FM,FS,FSS) 

C RESET D/A BOARD 
Y=2048 

STATUS=ALDV(0,Y) 

STATUS=ALTERM0 

WRITE(15,*)’Lambda = ’,LAMBDA, , R AND G = ’,R,G,’Payload = ’,PAYLOAD 
DO 200 M=0,NUM-1 
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200 WRITE(15,25)CURRENT(M+ 1),PSI(M),PSID(M+ 1),S(M + 1),PHI(M +1) 

25 FORMAT(5F10.4) 

CLOSE(15) 

END 

****************************** ******^****************************************** 
C Straight Sliding Mode Control 

SINCLUDE: ’ ATLDEFS. FOR ’ 

SINCLUDE: ’ ATLERRS.FOR’ 

PARAMETER (NUM = 1000) 

INTEGERS ADGAINS( 16),ADCHAN( 16),ICONFIG( 16),BASEADR,CDEVID, 
+CDEVF,LG,SCAN,DAVAL1,DA VAL2.Y,STATUS 
INTEGER*2 FH,FM,FS,FSS 
C 

REAL BHAT,FHAT,L,PI,K,R(3),IHAT 

REAL PSID(NUM),PSIVD(NUM),PSUERK(NUM),S(NUM),CURRENT(NUM) 

REAL PSI(0:NUM),PSIV(0:NUM),INTGS 

COMMON/CONFIG/ICONFIG 

EQUIVALENCE 

(ICONFIG(KCBASEADR) ,B ASEADR), (ICONFIG (KCDEVID), CDEVID) 

+ ,(ICONFIG(KCDEVFLAGS),CDEVFLG),(ICONFIG(KCSCAN),SCAN),(ICONFIG(KCC 
CHANNELS),CHAN) 

C 

DATA R/300.0,30.0,1.0/,SPEED,L/10.0,0.9985/ 

DATA SMALL,SMALFIL/2*0.0/,RATE,FREQ/2000.0,200.0/ 

DATA PSUERK/NUM*0.0/,PSIV(0),D20LD,VOLD,INTGS/4*0.0/ 

C UNCERTAINIES 

DATA BHAT,DB/1.1E4,295.0/,FHAT,DF/-5.7732E4,6.7114E4/ 

C 

PRINT *,’Please Enter Phi AND W (Real Number Please)!’ 

READ *,PHI,W 

PRINT *,’Please Enter Desired Trajectory "1" or "2" ’ 

READ VANS 
C DESIRED TRAJECTORY 
IF(IANS .EQ. 1) THEN 
DO 100 M = 1 ,NUM 
TIME=REA L(M)/FREQ 
IF (TIME .LE. 1.0) THEN 
PSID(M) =TIME 
PSIVD(M) = 1.0 

ELSE IF((TIME .GE. 1.0) .AND. (TIME .LE. 2.0)) THEN 
PSID(M) = 1.0 
PSIVD(M) =0.0 
END IF 
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IF((TIME .GT. 2.0) .AND. (TIME .LE. 2.4)) THEN 
PSID(M)=3.0-TIME 
PSIVD(M)=-1.0 
ELSE IF(TIME .GE. 2.4) THEN 
PSID(M)=0.6 
PSIVD(M)=0.0 
END IF 

100 CONTINUE 

ELSE IF (IANS .EQ. 2) THEN 
DO 500 M=1,NUM 
PSID(M) = 1.0 
PSIVD(M)=0.0 
PSUERK(M) =0.0 
500 CONTINUE 
END IF 
C 

M=0 

PI=4*ATAN(1.0) 

PAR=PI/1054.0 

R1B=1.0/(R(3)*BHAT) 

RlNOM=SQRT(R(l)**2+R(2)**2+R(3)**2) 

OPEN (UNIT = 15, FILE=’ C: \ ATL AB\TMP\CTRL2. M ’ .STATUS = ’NEW’) 

C INITIALIZE DATA ACQUSITION BOARD 
STATUS=ALINIT0 
STATUS=ALSB(1) 

STATUS=ALSF(RATE) 

STATUS=ALRSET() 

STATUS=ALGC(ICONFIG) 

STATUS=ALDV(0,Y) 

C 

WRITE(V)’ Are you ready (Ctrl-C) ? ’ 

WRITEC 1 ,*)’ When you are ready to go, press "1" and "return" ’ 

READ *,IANS1 
IF(IANS1 .EQ. 1) GOTO 7 
C BEGIN SIGNAL PROCESSING 
7 CALL GETTIM(IH,IM,IS,ISS) 

5 STATUS=ALA V(l.l.DAVALl) 

STATUS=ALAV(2,1 .DAVAL2) 

D1 =REAL(DAVALl-2048) 

D2=REAL(D AVAL2-2048) 

C FILTERED VERSION OF SMALL MOTION SIGNAL COMING OUT FROM STRAIN 
INDICATOR 

SMALFIL=0.94175 "-SMALFIL+0.029129 *(D20LD + D2) 

C CALIBRATED EQUATION BETWEEN DISPLACEMENT^m) AND SMALFIL 
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SMALL=-1.0*(SMALFIL+15.0) 

SMANG=SMALL/(100.0*L) 

C LARGE MOTION CALIBRATED EQUATION 
THETA=PAR*(D 1-30.0) 

C TOTAL ANGLE 

PSI(M) =THETA+SMANG 
C CONTROLLER DESIGN 
IF(M .EQ. 0) GOTO 1 
POSERR=PSI(M)-PSID(M) 

C FILTERED VERSION OF TOTAL VELOCITY, PSIV(M) 

VNEW=(PSI(M)-PSI(M-1))*FREQ 

PSIV (M) =0.975309*PSrV(M-1)+0.01234568 * (V OLD+VNEW) 
VELERR=PSI V (M)-PSIVD(M) 

S (M)=R(1) "‘POSERR+R(2)*VELERR 
INTGS=INTGS+(S(M)/FREQ) 

C NORMINAL VALUE 

IHAT=R( 1)"“(-VELERR)+R(3)*(PSUERK(M)-F T T. . / -((W , *2)" , INTGS) 
C TO FIND K 

K = (R1 NOM*(DF4-DB*(R 1 B*IH A. i))+SPEED)/(1.0-(DB*R IB)) 
LAMBDA = (SQRT( 1 + DB*R 1B)) * W 
PHI(M) = (K/ (2 *ZET A *LAMbD A)) * (DB “R1B) 

STD=S(M)/PHI(M) 

IF(ABS(S(M)) .LT. PHI(M)) SAT=STD 
IF(ABS(S(M)) .GE. PHI(M)) THEN 
IF(STD .GE. 0.0) SAT =1.0 
IF(STD .LT. 0.0) SAT=-1.0 
ENDIF 

CURRENT(M)=R1B"‘((IHAT - K*SAT)/1000.0) 

IF (CURRENT(M) .GE. 10.0E-03) CURRENT(M) = 10.0E-03 
IF (CURRENT(M) .LE. -10.0E-03) CURRENT(M)=-10.0E-03 
V OLT=CURRENT (M) "‘500.0 

C DIGITAL TO ANALOG (D/A) CALIBRATED EQUATION 
Y=NINT(((VOLT+9.9868)/4.8766)* 1000.0) 

STATUS=ALDV(0,Y) 

D20LD=D2 

VOLD=VNEW 

C 

1 M=M+1 

IF ( M .LE. NUM ) GOTO 5 
WRITE(*, *) * Game Over!’ 

CALL GETTIM(FH,FM,FS,FSS) 

C RESET D/A BOARD AND TERMINATION 
v =2048 

STATUS = ALDV(0,Y) 
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STATUS=ALTERMO 
WRITE(15,*)’%R = \R 
WRrrE(15,*)’%PHI= , ,PHI,’ W = ’,W 
DO 200 M=0,NUM-1 

200 WRITE(15,25)CURRENT(M + 1),PSI(M),PSID(M+1),PSIV(M),S(M+1) 

25 FORMAT(1X,5F12.4) 

CLOSE(15) 

END 

*************************************************************************** 


C Versatile Sliding Mode Control 
SINCLUDE: ’ ATLDEFS .FOR’ 

SINCLUDE: ’ ATLERRS. FOR’ 

PARAMETER (NUM=1000) 

INTEGER *2 ADGAINS(16),ADCHAN(16),ICONFIG(16),BASEADR,CDEVID, 

+CDEVF,LG,SCAN,DAVAL1 ,D A V AL2, Y, STATUS 
INTEGERS FH,FM,FS,FSS 
C 

REAL BHAT,FHAT,L,PI,K,R(3),G,IHAT 

REAL PSID(NUM), PSIVD(NUM), PSUERK(NUM), SDOT(NUM), CURRENT(NUM) 
REAL PSI(0:NUM),PSIV(0:NUM),INTGS,INTGERR 
COMMON/CONFIG/ICONFIG 
EQUIVALENCE 

(ICONFIG(KCBASEADR) ,B ASEADR), (ICONFIG(KCDEVID), CDEVID) 

+, (ICONFIG (KCDEVFLAGS), CDEVFLG), (ICONFIG(KCSC AN) ,SC AN), (ICONFIG(KCC 
CHANNELS),CHAN) 

C 

DATA R,G/300.0,30.0,1.0,1000.0/SPEED ,U 10.0,0.9985/ 

DATA SMALL,SMALFIL,INTGERR/3*0.0/,RATE,FREQ/2000.0,200.0/ 

DATA PSUERK/NUM*0.0/,PSIV(0),S,D20LD,VOLD,INTGS/5*0.0/ 

DATA SDOTO,SDOTN/2*0.0/ 

C UNCERTAINIES 

DATA BHAT, DB/1.1E4,295.0/, FHAT, DF/-5.7732E4,6.7114E4/ 

C 

PRINT *,’Please Enter Phi, Wn, AND Z (Real Number Please)!’ 

READ *,PHI,W,Z 

PRINT "‘/Please Enter Desired Trajectory "1" or "2" ’ 

READ "'.LANS 
C DESIRED TRAJECTORY 
IF(IANS .EQ. 1) THEN 
DO 100 M=1,NUM 
TIME=REAL(M)/FREQ 
IF (TIME .LE. 1.0) THEN 
PSID(M)=TIME 
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PSIVD(M) = 1.0 

ELSE IF((TIME .GE. 1.0) .AND. (TIME .LE. 2.0)) THEN 
PSID(M) = 1.0 
PSIVD(M)=0.0 
END IF 

IF((TIME .GT. 2.0) .AND. (TIME .LE. 2.4)) THEN 
PSID(M)=3.0-TIME 
PSIVD(M)=-1.0 
ELSE IF(TIME .GE. 2.4) THEN 
PSID(M) =0.6 
PSIVD(M)=0.0 
END IF 

100 CONTINUE 
ELSE 

DO 300 M = 1,NUM 
TIME=REAL(M)/FREQ 
PSID(M) = 1.0 
PSIVD(M)=0.0 
C PSIAD(M)=0.0 
PSUERK(M) =0.0 
300 CONTINUE 
ENDIF 
C 

M=0 

PI=4*ATAN(1.0) 

PAR=PI/1054.0 
R1B=1.0/(R(3)*BHAT) 

RlNOM=SQRT(R(l)**2-f-R(2)**2+R(3) , '‘*2) 

OPEN(UNIT=15,FILE=’C:\ATLAB\TMP\CTRL3.M’,ST ATUS =’NEW’) 
C INITIALIZE DATA ACQUSITION BOARD 
STATUS = ALINIT0 
STATUS=ALSB(1) 

STATUS = ALSF(RATE) 

STATUS = ALRSET0 
STATUS = ALGC(ICONFIG) 

STATUS=ALDV(0,Y) 

C 

WRITER,U’ Are you ready (Ctrl-C) ? ’ 

WRITE(*,*)’ When you are ready to go, press "1" and "return" ’ 

READ MANSI 
IF(IANS1 .EQ. 1) GOTO 7 
C BEGIN SIGNAL PROCESSING 
7 CALL GETTIM(IH,IM,IS,ISS) 

5 STATUS = ALAV(l,l,DAVAL1) 
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STATUS=ALAV(2,1 ,DAVAL2) 

D1 = REAL(D A V AL1-2048) 

D2=REAL(DAVAL2-2048) 

C FILTERED VERSION OF SMALL MOTION SIGNAL COMING OUT FROM STRAIN 
INDICATOR 

SMALFIL=0.94175*SMALFIL+0.029129*(D20LD+D2) 

C CALIBRATED EQUATION BETWEEN DISPLACEMENT(cm) AND SMALFIL 
SMALL=-1.0*(SMALFIL+15.0) 

SMANG=SMALL/(100.0*L) 

C LARGE MOTION CALIBRATED EQUATION 
THETA=PAR*(D 1-30.0) 

C TOTAL ANGLE 

PSI(M) =THETA+SMANG 
C CONTROLLER DESIGN (VERSATILE) 

IF(M .EQ. 0) GOTO 1 
POSERR=PSI(M)-PSID(M) 

INTGERR=INTGERR + POSERR/FREQ 
C FILTERED VERSION OF TOTAL VELOCITY, PSIV(M) 

VNEW=(PSI(M)-PSI(M- 1))*FREQ 

PSIV(M)=0.975309*PSIV(M-1)+0.01234568*(VOLD+VNEW) 

VELERR=PSIV (M)-PSIVD(M) 

SDOTN=R(l)*POSERR+R(2)*VELERR+G*INTGERR-Z*S 
SDOT(M)=0.94175*SDOT(M-1)+0.029129*(SDOTO+SDOTN) 
S=S+SDOT(M)/FREQ 
C NORMINAL VALUE 

IHAT=R(l)*(-VELERR)+R(3)*(PSUERK(M)-FHAT)-G*POSERR+Z*SDOT(M) 

+-((W**2)*S) 

CTO FIND K 

K=(R1NOM*(DF+DB*(R1B*IHAT))+SPEED)/(1.0-(DB*R1B)) 

LAMBDA=(SQRT(1 +DB*R1B))*W 

PHI(M)=(K*(DB*R1B))/((2*ZETA*LAMBDA)+(DB*R1B*Z)) 

C 

STD=SDOT (M)/PHI(M) 

IF(ABS(SDOT(M)) XT. PHI(M)) SAT=STD 
DF(ABS(SDOT(M)) .GE. PHI(M)) THEN 
EF(STD .GE. 0.0) SAT=1.0 
IF(STD .LT. 0.0) SAT=-1.0 
ENDIF 

CURRENT(M)=R1B*((IHAT - K*SAT)/1000.0) 

IF (CURRENT(M) .GE. 10.0E-03) CURRENT(M) = 10.0E-03 
IF (CURRENT(M) .LE. -10.0E-03) CURRENT(M)=-10.0E-03 
VOLT=CURRENT(M) *500.0 

C DIGITAL TO ANALOG (D/A) CALIBRATED EQUATION 
Y=NINT(((VOLT+9.9868)/4.8766)* 1000.0) 


114 









STATUS=ALDV(0,Y) 

D20LD=D2 

VOLD=VNEW 

SDOTO=SDOTN 

C 

1 M=M + 1 

IF ( M .LE. NUM ) GOTO 5 
WRITER, *) ’Game Over !’ 

CALL GETTIM(FH,FM,FS,FSS) 

C RESET D/A BOARD AND TERMINATION 
Y=2048 

STATUS=ALDV(0,Y) 

STATUS=ALTERM0 
WRITE(15,*)’%R=\R 
WRITE(15,*)’%G = \G 

WRITERS,*)’% PHI = ’,PHI,’ W= ’,W,’ Z = ’,Z 
DO 200 M=0,NUM-1 

200 WRITE(15,25)CURRENT(M+l),PSI(M),PSID(M + l),PSIV(M),SDOT(M+l) 
25 FORMAT(lX,5F12.4) 

CLOSE(15) 

END 
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